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PREFACE 


Verification  of  a  coiq>rehen8lve  test  ban  treaty  involves  certain 
technical  issues.  This  report,  sponsored  by  the  U.S.  Ams  Control 
and  Dlsamanent  Agency,  is  part  of  a  continuing  research  program  on 
such  Issues.  It  presents  the  analysis  and  code  documentation  for 
SNAP/D,  a  computer  model  for  assessing  the  ability  of  a  network  to 
detect  and  locate  seismic  events.  The  module  in  SNAP/D  that  calcu¬ 
lates  the  uncertainty  in  network  location  estimates  was  developed  by 
Systems,  Science  and  Software  under  sponsorship  of  the  Defense  Advanced 
Research  Projects  Agency. 
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I.  INTRODUCTION  AMD  SUMMARY 


This  report  is  a  user's  manual  for  SNAP/D  (Seismic  Network 
Assessment  Program  for  Detection),  a  computer  model  designed  ^  ~ 

- - BaaaayplTi  ftn  (Pgl>:>- f ny  aSSeSSlng  the  ability 

of  a  network  to  detect  and  locate  seismic  events  (l.e.,  earthquakes 
or  explosions).  Part  I  presents  the  analytical  structure  of  the 

X _ _ 

model,  whereas  Part  II  documents  the  computer  code.  V 

PSR  plans  to  develop  another  model,  SNAP/I,  for  assessing  the 
ability  of  a  seismic  network  to  identify  earthquakes  and  explosions 
[Clervo,  1978].  S^^AP/D  will  be  used  as  a  starting  point  for  SNAP/I, 
inasmuch  as  a  seismic  event  must  be  detected  before  it  can  be  iden¬ 
tified.  Therefore,  SMAP/D  includes  many  features  that  will  be  used 
by  SNAP/I  to  distinguish  earthquakes  from  explosions.  For  example, 
the  user  can  specify,  for  any  SNAP/D  run,  either  an  explosion,  a 
shallow  earthquake,  or  a  deep  earthquake— each  of  which  propagates 
waves  at  relatively  different  average  amplitudes. 

In  addition  to  its  eventual  role  in  SNAP/I,  SNAP/D  is  expected 
to  be  of  immediate  use  to  analysts  concerned  with  test  ban  treaty 
verification.  It  is  intended  to  replace  NETWORTH  [Wirth,  1977],  the 
computer  model  that  has  gained  wide  acceptance  in  the  verification 
community  for  assessing  network  performance,  SNAP/D  has  a  much 
broader  scope  than  NETHDRTH.  For  example,  SNAP/D 

•  "Propagates"  (models  the  propagation  of)  up  to  10  seismic 
waves  in  a  single  run. 

•  Calculates  wave  attenuation  and  travel  time  as  functions  of 
focal  depth,  regional  media  characteristics,  and  event  type. 

•  Allows  user  specification  of  multiwave  network  detection 
criteria. 

•  Accosnodates  advanced  selsmological  techniques  for  deter¬ 
mining  the  uncertainty  in  estimates  of  focal  depth  and  epi¬ 
center  location,  such  as 
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•  Master-event  calibration  of  P  travel  tines, 

•  S-P  origin-time  fix, 

•  pP  depth  estimates  with  step-out, 

•  Lg  back  azimuth  estimates. 

•  Calculates  the  relative  Importance  of  individual  waves  and 
stations  to  hypocenter  estimates. 

•  Easily  Incorporates  new  waves  or  techniques  because  of  its 
modular  design'. 

•  Allows  the  user  to  initiate  runs  in  either  an  interactive 
or  a  batch  mode. 


The  seismic  waves  (or  phases)  propagated  in  the  model  include  P, 
Pgt  Lg,  Sn,  S,  pP,  and  Rayleigh,  with  provision  for  propagating 
four  additional  phases.  The  existence  and  attenuation  of  each  phase 
with  eplcentral  distance  is  keyed  to  event  type  (earthquake  or  explo¬ 
sion),  focal  depth,  and  event  size  (generally  given  in  m^  units). 

In  addition,  regional  media  characteristics  along  the  propagation 
path  are  used  to  compute  the  expected  amplitude  of  each  phase  observed 
at  a  given  station. 

At  the  start  of  each  run,  the  user  selects  either  an  explosion 
at  zero  depth,  a  shallow  earthquake  at  15  km  depth,  or  a  deep  earth¬ 
quake  at  100  km  depth.  The  selected  event  is  assumed  to  occur  over 
one  or  several  grids  of  hypothetical  epicenters.  In  its  most  common 
operating  mode,  SMAP/D  computes  the  minimum  magnitude  in  m^  units,  at 
each  epicenter,  that  provides  a  given  probability  that  the  selected 
network  of  stations  will  meet  the  selected  multiwave  detection  cri¬ 
teria,  For  example,  at  a  given  epicenter,  SNAF/D  could  compute,  with 
0.90  probability,  the  smallest  explosion  (in  units)  that  would 
result  in  at  least  a  four-station  F-wave  observation  or  a  two-station 


F-  and  Lg-wave  observation. 

SNAF/D  can  also  calculate,  for  each  hypothetical  epicenter,  the 
eplcentral  error  ellipse  and  depth  confidence  Interval  at  any  user- 
specified  level.  The  generalized  inversion  technique  used  to  compute 
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those  statistical  parameters  Incorporates  multiwave  arrivals  at  all 
stations  as  well  as  specialized  propagation-path  calibration  proce¬ 
dures  specified  by  the  user.  The  procedures  include  master-event 
calibration  of  F-wave  travel  times,  pP  depth  estimates  with  step-out, 
and  S-P  origin-time  fix.  Back  azimuth  estimates  from  Lg-wave  obser¬ 
vations  are  also  included  in  the  error  calculations.  The  relative 
"importance"  of  individual  waves  and  stations  to  the  precision  of 
hypocenter  estimates  at  any  epicenter,  or  set  of  epicenters,  may  also 
be  computed. 
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PART  I 

AMALYTICAL  MODEL 
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II .  MACROSTRUCTURE 


This  section  specifies  the  waves  propagated  In  SNAP/D,  describes 
the  SNAP/D  modules,  and  presents  a  flowchart  of  the  data  passed  between 
modules.  The  modular  design  of  SNAP/D  has  proved  valuable  In  concep¬ 
tualizing  and  describing  model  operations.  However,  some  modules  con¬ 
tain  several  subroutines,  and  some  subroutines  contain  more  than  one 
module.  Part  I  presents  SNAP/D  from  a  modular  point  of  view,  whereas 
Part  II  describes  the  operation  of  the  subroutines. 

I 

SEISMIC  WAVES 

Table  1  lists  the  waves  propagated  In  SMAP/D  along  with  their 
primary  functions  (detection,  location,  or  Identification) ,  regional 
propagation  path  corrections,  and  seismic  source.  The  P  and  Pn  as  well 
as  the  S  and  Sn  phases  have  been  combined  Into  one  wave  for  shallow 
events  (explosion  or  shallow  earthquake),  because  of  nearly  nonover¬ 
lapping  ranges  and  similar  functions.  The  reflected  phase  pP  repre¬ 
sents  both  pP  and  sF  In  SNAP/D,  since  both  serve  a  similar  function 
In  depth  estimation  but  are  rarely  observed  by  the  same  station  for  a 
given  event.  Furthermore,  because  sP  Is  more  sensitive  than  pP  to 
variations  In  depth,  we  can  conservatively  assess  the  role  of  the  re¬ 
flected  phases  In  depth  estimation  by  using  the  pP  travel-time  table 
to  represent  both  pP  and  sP  travel  times. 

Although  any  wave  listed  in  Table  1  may  be  Included  In  the  net¬ 
work  detection  criteria,  detection  Is  listed  as  a  primary  function 
only  for  waves  that  have  either  been  traditionally  used  for  network 
detection  or  are  proposed  for  such  use.  Since  network  detection  cri¬ 
teria  have  been  based  on  the  ability  of  the  network  to  locate  an  event, 
any  wave  In  SNAP/D  that  Is  listed  for  detection  is  also  listed  for 
location  (i.e.,  used  to  calculate  the  precision  of  epicenter  location 
and  focal  depth  estimates).  Further,  since  we  plan  to  Incorporate 
SNAP/D  into  an  Identification  model  (SNAP/I),  we  have  made  the 
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provision  for  propagating  waves  that  will  be  used  for  identifying 
seismic  events.  The  Rayleigh  wave  used  in  the  discriminant  is 

one  such  wave,  but  SNAP/D  can  also  propagate  four  additional  phases 
for  use  in  seismic  discrimination. 

In  addition  to  the  traditional  least  squares  "trlangularization" 
process  based  on  arrival  times.  Table  1  indicates  the  advanced  loca¬ 
tion  techniques  modeled  in  SMAP/D.  Thus,  for  a  user-specified  set  of 
stations  calibrated  with  respect  to  any  epicenter,  SNAP/D  accounts 
for  both  master-event  (i.e.,  large  earthquake)  calibration  of  P 
arrivals  and  Vp/v^  path  calibration  for  use  in  determining  origin  time 
from  the  difference  between  S  and  P  arrivals.  The  sets  of  calibrated 
stations  for  any  epicenter  are  specified  separately  for  each  type  of 
calibration  (master-event  or  S-P),  and  each  source  depth  (shallow  or 
deep).  The  Lg  phase  provides  back-azimuth  rather  than  arrival-time  in¬ 
formation  for  calculating  the  precision  of  hypocenter  estimates.  SNAP/D 
also  models  the  "step-out"  procedure  for  use  with  the  reflected  phase  pP 
for  depth  estimation. 

At  the  start  of  each  run,  the  user  can  choose  from  three  event- 
type/  focal-depth  comblnatlond:  (1)  an  explosion  (X)  at  zero  depth; 

(2)  a  shallow  earthquake  ((^)  at  15  km  depth;  or  (3)  a  deep  earth¬ 
quake  (Qq)  at  100  km  depth.  Table  1  summarizes  the  general  effect  of 
that  choice  on  the  waves  propagated  in  the  model.  As  shown  in  the 
table,  SNAP/D  usually  performs  the  amplitude  corrections  for  differ¬ 
ent  events  by  adding  factors  to  magnitudes  computed  from  standard 
attenuation  tables.  However,  a  new  travel-time  table  is  required  for 
P-wave  propagation  from  a  deep  earthquake. 

In  general,  the  probability  that  a  given  station  will  observe  a 
given  wave  from  any  epicenter  is  computed  from  the  expected  amplitude 
and  standard  deviation  of  the  wave  and  associated  noise  at  that  sta¬ 
tion.  The  expected  amplitude  of  the  wave  is  found  from  a  standard 
attenuation  table  (or  tables)  corrected  for  event  size  (in  mj^  units), 
event  type,  and  regional  media  characteristics  along  the  propagation 
path.  However,  propagation  path  corrections  are  made  only  when  the 
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source  is  shallow  (X  or  Qg)  and  the  distance  from  epicenter  to  station 
Is  less  than  20  deg.  Propagation  path  corrections  begin  by  computing 
the  relative  portion  of  the  path  that  lies  In  tectonic  or  Intermediate 
(as  opposed  to  stable)  rock  from  a  user-specified  regional  media  grid. 
That  computation  Is  then  used  to  obtain  the  expected  wave  amplitude 
for  a  given  eplcentral  distance  by  Interpolating  between  stable  and 
tectonic  attenuation  t  '  les.  Similarly,  SI  .P/D  accounts  for  Lg  phase 
attenuation  by  large  bodies  of  water  and  mountain  ranges  by  means  of 
a  user-specified  "Lg  squelch"  grid. 

INFORMATION  ROUTING 

Figure  1  shows  the  flowchart  of  SNAP/D.  The  control  Inputs  serve 
two  major  purposes:  first,  they  define  the  nature  of  the  run  (detec¬ 
tion  or  location,  explosion  or  earthquake,  focal  depth,  etc.);  second, 
they  preset  all  other  Inputs  and  modules  so  that  unnecessary  functions 
are  not  performed  during  a  model  run.  Detection  inputs  define  the 
network  detection  criteria  and  such  associated  parameters  as  the  proba¬ 
bility  threshold  for  output  magnitude  contours.  Epicenter  Inputs 
establish  the  grid  of  test  epicenters.  Regional  inputs  provide  the 
grid  of  regional  corrections  for  local  wave  propagation.  Station 
Inputs  Include  station  location,  reliability,  and  noise  statistics. 

The  calibration  and  importance  inputs,  respectively,  specify  the 
epicenter/station  pairs  for  travel-time  calibration  and  Importance- 
matrix  calculations.  The  wave  Inputs  primarily  consist  of  amplitude 
and  travel-time  tables  and  associated  statistics.  A  detailed  descrip¬ 
tion  of  all  Inputs  Is  presented  In  the  next  section. 

The  calculation  modules  shown  In  Fig.  1  have  been  chosen  for 
easy  conceptualization  of  SNAP/D.  Each  module  may  contain  several 
subroutines  of  the  computer  program  or  vice  versa  (as  discussed  in 
Part  II  of  this  report).  Each  module  function  Is  sumnarlzed  below 
and  detailed  In  Sec.  IV, 


•  TRANS— Translates  the  multiwave  detection  criteria  Into  two 

types  of  canonical  calculations,  then  uses  the  calculations 
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to  compute  the  network  detection  probability.  TRANS  also 
initiates  printing  of  an  error  message  if  the  user  enters 
inadmissible  detection  criteria.  Network  detection  criteria 
are  entered  in  a  simple,  intuitive  format. 

•  NPROB-— Performs  the  canonical  probability  calculations  re¬ 
quired  by  TRANS. 

•  EGRID — Uses  the  epicenter  inputs  to  set  up  the  grid  of  test 
epicenters,  and  increments  the  program  from  epicenter  to 
epicenter , 

•  BINMB— Increments  the  program,  through  magnitude  values,  in 
a  binary  search  for  the  magnitude  that  results  in  network 
detection  at  the  user-entered  probability  threshold.  BINMB 
performs  these  operations  for  each  test  epicenter. 

•  GEQM — Calculates  the  great-circle  geometric  parameters  re¬ 
quired  by  the  model.  GEOM  also  computes  the  regional  weight 
from  epicenter  to  station  which  determines  the  relative  por¬ 
tions  of  the  propagation  path  in  stable  or  tectonic  crustal 
media. 

•  WAMF — Normalizes  all  waves  propagated  in  the  model  to  the 
size  of  the  event  (in  m^^  units),  and  corrects  them  for  event 
type  (explosion  or  earthquake)  and  depth  of  focus.  WAMP 
performs  this  function  by  calculating  the  expected  log  ampli¬ 
tude  for  each  wave  observation. 

•  SPROB — Calculates  probability  that  station  1  will 

observe  wave  k  propagated  from  an  event  at  epicenter  J,  for 
all  (Ijk)  relevant  to  the  particular  jua. 

m  ATSD — Computes  the  arrival-time  standard  deviation  for  each 

wave.  The  computation  accounts  for  the  advanced  techniques 
used  to  obtain  location  estimates  such  as  master-event  cal¬ 
ibration  of  P-wave  arrivals  and  S-P  time  fix. 

•  SUBSET — Calculates  all  probabilities  associated  with  station 
subsets.  Thus,  SUBSET  calculates  the  probability  that  the 
network  has  met  the  pP  step-out  criteria.  For  SNAP/I,  SUBSET 
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vould  calculate  the  probability  that  the  network  has  satis¬ 
fied  the  first-notion  discriminant  criteria. 

•  LOCATE — Calculates  the  area  of  the  epicentral  estimate  error 
ellipse,  the  length  of  the  depth  confidence  interval,  and 
the  Importance  of  vave/station  combinations  at  user-specified 
epicenters . 

•  PRINT— Controls  and  formats  the  output  of  SNAP/D. 

SNAP/D  performs  two  largely  independent  sets  of  computations. 

If  the  user  selects  a  detection  contour  printout  alone,  then  the  ATSD 
and  LOCATE  modules  are  inoperative  for  that  run.  On  the  other  hand, 
if  the  user  selects  only  a  location  contour  printout  for  a  specified 
event  size,  then  TRANS,  NPROB,  and  BINMB  are  inoperative  for  that  run. 
If  the  user  selects  both  sets  of  computations,  he  may  perform  loca¬ 
tion  calculations  for  each  epicenter  at  the  detection  threshold  mag¬ 
nitude  (in  units) . 
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III,  INPUTS 


This  section  sunmarlzes  the  Inputs  to  SNAF/D.  Part  II  details  the 
interactive  format  in  which  they  are  used.  However,  the  user  should 
fsailiarize  himself  with  the  material  presented  here  before  attempting 
to  run  SNAP/D,  since  the  interactive  prompts  in  Part  II  are  too  brief 
to  completely  describe  the  function,  units,  range,  etc.,  of  each  input. 

The  input  categories  and  file  structure  described  below  conform 
to  the  convention  used  for  interactive  prompting.  When  a  file  is  in¬ 
dicated,  the  corresponding  data  may  be  stored  and  recalled  for  use 
with  or  without  review  or  updating.  The  input  categories  are 


•  Control  inputs — Nonnumerlcal  information  that  defines  the 
nature  of  the  run. 

•  Detection  input  file— ’A  short  file  that  provides  the  multi¬ 
wave  network  detection  criteria  and  associated  parameters 
for  runs  computing  detection  contours. 

m  Topographic  inputs— Inputs  that  require  specification  on  a 
world  map.  Topographic  inputs  are  filed  in  the  following 
five  categories: 

•  Station  file;  Locates  and  statistically  describes  each 
station  in  the  seismic  network. 

•  Epicenter  grid  file:  Provides  the  set  of  test  epicenters, 
which  may  or  may  not  be  broken  into  discontiguous  blocks. 

•  Eplcenter/statioh  calibration  file:  Lists  the  set  of  sta¬ 
tion  indices  that  are  master-event  and  S-P  calibrated 
with  respect  to  each  epicenter  index.  For  aonsiatent  in¬ 
dexing,  each  epioenter/atation  oalibration  file  muat  be 
aroaa-referenoed  to  a  given  pair  of  atation  and  epicenter 
grid  filea, 

•  Regional  grid  file:  Provides  stable/ tectonic  and  Lg 
squelch  weights  with  a  simplified  default  mechanism. 

■  Importance-matrix  file:  Determines  the  sets  of  epicenters 
for  which  wave/station  importance  values  will  be  calculated. 


14 


AC90C114 


•  Wave  Inputs — Inputs  that  pertain  to  the  waves  listed  in 
Table  1  (p.  8).  The  inputs  for  each  wave  are  filed  under 
the  name  of  the  wave.  Each  wave  input  file  contains  sub¬ 
files  for  each  attenuation,  travel-time,  or  travel-time 
standard-deviation  table,  as  required. 

Recall  that  each  SNAP/D  run,  as  established  by  the  control  in¬ 
puts,  is  for  a  specific  event/depth  combination:  an  explosion  (X) 
at  0  km  depth,  a  shallow  earthquake  (Qg)  at  15  km  depth,  or  a  deep 
earthquake  (Q^)  at  100  km  depth.  Table  1  (p.  8)  shows  that  some 
waves  are  not  generated  by  a  particular  event.  Thus,  certain  inputs 
are  required  for  certain  event  types  and  not  others — as  indicated  in 
the  descriptions  provided  below.  Furthermore,  input  values  may  de¬ 
pend  on  event  type;  thus,  the  user  must  indicate  the  relevant  event 
in  each  file  description. 

In  presenting  the  inputs  below,  we  have  followed  the  same  se¬ 
quence  and  terminology  used  for  the  SNAP/0  interactive  prompting. 
However,  the  discussions  have  been  expanded  Co  give  Che  user  a  better 
understanding  of  the  nature  and  function  of  each  input. 

CONTROL  INPUTS 

Select  one  option  for  each  of  the  following  inputs: 

1.  Calibration  mode: 

D,  detection  only 

L,  location  only 

DL,  detection  and  location 

2.  Event  type: 

X,  explosion  at  0  km 

Qg,  shallow  earthquake  at  15  km 

Qq,  deep  earthquake  at  100  km 
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3.  Detection  output  option  for  D  and  DL  runs: 

MAG,  magnitude  contours  for  fixed  network-detection 
probability 

PR,  probability  contours  for  fixed  magnitude  value 

4.  will  importance-matrix  calculations  be  required? 

Yes  or  no. 

DETECTION  INPUT  PILE 

Since  the  detection  Inputs  are  filed,  the  user  may  create  new 
files  or  review  and  update  old  ones.  Detection  inputs  are  required 
for  only  D  or  DL  runs. 

1.  Name  and  description  of  the  file. 

2.  The  detection  criteria.  (User  is  alerted  at  this  point  if 
an  inadmissible  criterion  is  entered.)  For  exaiaple  (in 
interactive  format): 

a.  ((P  +  PG)*LG)/2  +  (P)/4,  with  the  literal  meaning,  (P 
or  Pg)  and  Lg  observed  by  at  least  two  stations,  or  P 
observed  by  at  least  four  stations. 

b.  (P)/4,  meaning,  at  least  a  four-station  P-wave  obser¬ 
vation. 

c.  (R)/2,  meaning,  at  least  a  two-station  Rayleigh-wave 
observation. 

If  output  option  MAG  was  ohosen  in  control  inputs^  the  user 
nou  enters 

3.  P^,  the  probability  threshold  for  magnitude  contours  in 
either  units  for  multiwave  or  P-wave-only  detection  cri¬ 
teria  (examples  a  and  b  above^,  or  in  m^  units  for  single¬ 
wave  detection  criteria  involving  wave  k  (e.g.,  magnitude 
contours  in  units,  such  as  example  c  above). 

4.  the  lowest  magnitude  used  in  the  binary  search  in  m^^ 
or  mj^  units,  as  above. 

5.  the  highest  magnitude  used  in  the  binary  search  in 
or  m^^  units,  as  above. 
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If  output  option  PR  waa  ohoaen  in  aontrol  inputs ^  the  user 
noa  enteva 


6.  the  magnitude  threshold  in  either  m^  units  (examples  a 

and  b  above),  or  units  (exaiq>le  c  above). 


TOPOGRAPHIC  IMPUTS 

Topographic  inputs  consist  of  one  file  from  each  of  the  follow¬ 
ing  five  categories:  station  files,  epicenter  grid  files,  epicenter/ 
station  calibration  files,  regional  grid  files,  and  importance-matrix 
files.  Users  should  have  a  map  of  the  relevant  area  available  before 
creating  any  of  those  files,  since  correspondence  of  station,  epicenter, 
and  regional  grid  indexing  with  geographical  coordinates  is  essential 
for  ;  roper  operation  of  SNAP/D. 


^tatiofc  Pile 

The  first  seven  of  the  following  inputs  apply  to  all  stations  in 


tha  network.  Amplitude  units  are  0-P  m|i. 


1. 

2. 

3. 

4. 

5. 

6. 

7. 


Name  and  description  of  the  file. 

N,  the  number  of  stations  in  the  seismic  network. 

the  mean  of  the  lognormally  distributed  noise  level  for 
wave  k  at  all  stations,  k  ■  1,  ...,  10. 

a^,  the  standard  deviation  (s.d.)  of  noise  amplitude 

for  wave  k  at  all  stations. 

the  s.d.  of  signal  amplitude  for  wave  k  at  all 

stations  for  a  given  magnitude  nij^,  in  the  magnitude  units 
for  wave  k. 


Obk>  the  additional  s.d.  of  IoRj^q  signal  amplitude  for  wave  k 
at  all  stations,  given  an  value.  This  input  is  not  re¬ 
quested  for  the  P  wave  (k  ■  1). 


tj^,  the  slgnal-to-noise  ratio  (SNR)  required  to  observe  wave  k 
at  all  stations. 
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The  Inputs  for  station  1,  1  *  1,  ,  N,  are 

8.  Naae  and  station,  four  letters  maximum. 

9.  (6^,  the  latitude  and  longitude  of  station  1.  6^  Is 

In  degrees  north  or  south  of  the  equator,  and  Is  In 
degrees  east  or  west  of  Greenwich  meridian. 

10.  R^,  the  reliability  of  station  1.  This  Is  the  probability 
that  station  1  Is  operational  at  any  time. 

11.  the  mean  noise  level  for  wave  k  at  station  1. 

If  ■  0,  then  noise  level  is 

if  Wniij  1*  Oj  then  noise  level  is 

12.  ^°^10  wave  k  at  station  1. 

If  ■  0,  then  noise  s.d,  is 

if  <Jnik  ^  then  noise  s,d.  is 

13.  ^^^10  wave  k  at  station  1. 

If  ensile  ■  0,  then  signal  s,d.  is  Ogiy* 
if  Ogj^jj  ^  0,  then  signal  s,d,  is 

14.  r^^,  the  SNR  required  to  observe  wave  k  at  station  1. 

If  rik  ■  0,  then  SNE  is  r^; 
if  rj^jj  ^  0,  then  SNE  is  rj^|j. 

15.  s^^,  the  station-dependent  arrival-time  s.d.  for  wave  k  at 
station  1,  k  •  1,  2,  5. 

16.  the  s.d.  of  the  Lg  (k  ■  3)  back  azimuth  estimate  error 
at  station  1. 

Epicenter  Grid  File 

1.  Name  and  description  of  the  file. 

2.  H,  the  number  of  epicenter  blocks.  Blocks  may  be  used  to 
provide  discontiguous  epicenter  grids  In  a  single  run  (e.g.. 
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NTS  and  USSR) ,  or  to  provide  additional  density  of  test 
epicenters  within  a  large  geographical  region, 

3.  the  number  of  epicenter  latitudes  for  block  h. 

4.  the  number  of  epicenter  longitudes  for  block  h. 

•  is  the  total  number  of  epicenters  in 
block  h. 

H 

•  J  ”  cannot  exceed  400. 

h-1 

5.  the  starting  latitude. 

6.  the  latitudinal  distance  between  epicenters. 

7.  the  starting  longitude. 

8.  the  longitudinal  distance  between  epicenters. 
Epicenter  indexing  proceeds  as 


ie<».  *  A*<») . 

(e‘».  ♦<»  .  [e(»  *  A6'», 

X  1  9  X 


9j^  J »  •  ♦  •  > 


.(1) 


(2) 


1, 


Epicenter/Station  Calibration  File 

1.  Name  and  description  of  the  file.  To  avoid  indexing  errors, 
the  description  should  cross-reference  the  station  and  epi¬ 
center  grid  files  for  the  relevant  calibration  Inputs. 

2.  tbe  magnitude  correction  for  station  1,  epicenter  J, 

and  wave  k.  Defaults  to  -  0  for  nonentries. 

3.  Qj,  the  set  of  stations  (i.e.,  station  Indices)  that  are 
master-event  calibrated  with  respect  to  an  event  at  epicen¬ 
ter  j. 

4.  i|>j,  the  set  of  stations  that  are  S-P  calibrated  with  respect 
to  an  event  at  epicenter  J. 
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Since  and  the  seta  of  calibrated  etationa  far  epicenter 
j,  may  depend  on  the  focal  depth  of  the  events  the  user  must 
indicate  in  the  calibration  file  name  or  description  whether 
the  file  pertains  to  a  shallow  ^X  or  Q^J  or  deep  ^QJ^J  event. 


Regional  Grid  File 

Each  regional  grid  file  covers  the  entire  globe.  Thus,  the  ssme 
file  can  be  used  %rlth  any  epicenter  file.  However,  the  user  Is  asked 
to  enter  weights  for  only  nonstable  (v^  ^  0)  and  Lg  squelching 
(Vj^  ^  0)  grids.  All  other  grid  weights  default  to  w^  *  0  and  v^  >  0. 


1.  Name  and  description  of  the  file. 

2.  the  grid  dimension. 

•  Individual  grids  have  equiangular  sides  so  that 

A6  ■  •  Aj^. 

•  In  general,  2*  <  Aj^  a:  10* . 

3.  Ad^,  the  size  of  polar  caps  in  degrees  latitude. 

•  Because  of  the  increasing  density  of  grids  near  the 
poles,  the  polar  caps  are  assigned  a  single  regional 
weight . 

•  The  starting  latitude  for  generating  the  first  regional 

grid  is  6j(2  ~  ~  B^’uth  of  the  equator.  The  start¬ 

ing  longitude  is  *  180*. 

•  180*  -  2A6^  and  360*  must  be  multiples  of 

4.  w„,  the  regional  weight  for  the  north  polar  cap,  0  $  w„  ^  1, 

£1  M 

5.  Wg,  the  regional  weight  for  the  south  polar  cap,  0  jC  w^  <  1. 

6.  ^  *  1,  •  •  •  t 

grid  excluding  the  polar  caps. 


the  weight  assigned  to  each  regional 


■  0  s  w|l  <  1,  where  w^  •  0  implies  a  purely  stable  geology 

for  the  %th  grid,  and  W£  ■  1,  purely  tectonic.  If  no 
weight  is  entered,  the  program  defaults  to  W£  ■  0. 

•  Ljj  -  (360*/Aj^)(180*  -  2Ae^/Aj). 
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7 •  •••» 

•  >  0  or  1,  where  V£  >  1  implies  that  Lg  waves  are 

squelched  in  the  Ith  grid,  and  V£  >  0  iisplies  no 
squelching.  The  interactive  prompts  for  SNAP/D  simply 
request  the  Indices  of  the  grid  for  which  V£  >  1. 


L^,  the  squelch  weight  for  Lg  waves. 


Importance^latrix  File 

If  importance-matrix  calculations  were  requested,  the  user  must 
now  specify  the  epicenters  to  which  the  calculations  will  apply.  An 
output  matrix  is  provided  for  each  input  set  of  epicenters. 


1.  Name  and  description  of  the  file. 

2.  r  ,  m  -  1,  ...,  M,  the  set  of  epicenter  indices  for  the 

IB 

mth  output  Importance  matrix. 

•  As  a  special  case,  the  user  may  request  an  importance- 
matrix  output  for  all  epicenters  without  actually  enum¬ 
erating  them. 


WAVE  INPUTS 

The  wave  inputs  consist  of  a  common-wave  file  and  a  file  for 
each  individual  wave.  The  common-wave  file  requires  the  following 
entries : 


1.  Name  and  description  of  the  file. 


2.  value  used  for  all  location  and  Importance 

calculations. 


If  mjjL  ■  0,  t?ien  SNAP/D  defaults  to  using  threshold  detection 
magmtude  for  calculating  wave  observation  probabilities  at  each 
station  for  a  given  epicenter.  However,  it  is  inadmissible  to 
set  m^jj^  •  0  for  an  L  run,  or  for  D  or  DL  nets  not  using  the 
P  wave  in  the  detection  criteria.  Also,  m|,L  should  not  be  set 
to  zero  if  a  fixed-magnitude  run  is  requested  (i.e.,  the  output 
option  PR  was  chosen  in  control  inputs). 


3.  P^,  the  confidence  probability  associated  with  error  ellipse 

and  depth-interval  calculations,  in  multiples  of  5  percent. 
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4.  The  Indices  of  the  waves  that  are  to  be  excluded  from  loca¬ 
tion  and  importance  calculations  (the  effect  of  excluding  a 
wave  is  simply  to  zero  out  its  amplitude  in  SNAP/D) . 

5.  The  names  of  the  wave  files  for  waves  that  will  be  used  for 
either  location  or  detection. 

The  individual  wave  file  Inputs  are  given  below.  All  attenuation, 
travel-time,  and  travel-time  s.d.  tables  are  in  files  that  have  en¬ 
tries  in  the  standard  (l.e.,  NETWORTH  [Vlrth,  1977])  format. 


P-Wave  Fj ’ * 

1.  *]ame  and  description  of  the  file. 

2.  Stable  media  attenuation  file  with  table  entries  starting 

at  A  -  0*. 

•  The  P-wave  stable  media  attenuation  file  name  should 
indicate  whether  the  attenuation  table  applies  to  a 
shallow  or  deep  event,  since  the  shallow  attenuation 
table  starts  with  the  Pn  phase. 

3.  Tectonic  media  atten\iation  file  for  A  £  20”. 

4.  Travel-time  file  for  starting  at  A  ■  0®. 

•  The  travel-time  file  name  should  indicate  whether  the 
table  is  applicable  to  a  shallow  or  deep  event. 

5.  Travel-time  s.d.  file. 

6.  k  *  1,  P-wave  arrival-time  reading  error  s.d. 

7.  Vp(z),  the  P-wave  velocity  at  z  ■  0,  15  or  100  km. 

•  Since  different  travel-time  tables  apply  to  shallow  and 
deep  events,  and  vp  is  a  function  of  depth,  the  user 
must  indicate  in  the  P-wave  file  name  or  description 
whether  the  file  applies  to  an  X,  Qg,  or  Qq  event.  For 
the  Qq  event,  the  tectonic  media  attenuation  file  is  not 
required . 


Pg-Wave  File 

The  Pg  wave  is  not  generated  by  a  event. 
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1.  Nane  and  daacrlptlon  of  the  file, 

2.  Stable  media  attenuation  file. 

3.  Tectonic  media  attenuation  file  (A  £  20”). 

4.  Travel-time  file. 

5.  Travel-time  a.d.  file, 

6.  event/depth  normalization  constants  for  stable 
(a  *  S)  and  tectonic  (a  >  T)  med:'a. 

7.  magnitude  normalization  constants  for  stable 
(a  "  S)  and  tectonic  (a  >  T)  media. 

*  The  use  of  [KE]^a2  explained  In  Sec.  IV, 

during  discussion  of  the  UAMP  module. 


The  Lg  wave  Is  not  generated  by  a  Q.  event.  Since  the  Lg  event 

tn\  ^ 

factor  (KE]^“2  Is  different  for  the  X  and  events,  the  Lg-wave  file 
name  or  description  should  Indicate  which  of  these  events  is  desired. 

1.  Name  and  description  of  the  file. 

2.  Stable  media  attenuation  file. 

3.  Tectonic  media  atteniiatlon  (A  <  20*) 

4.  a  -  S  or  T. 

5.  [KMlJSj,  a  -  S  or  T. 

S-Wave  File 

The  S  wave  Is  generated  by  all  events,  but  Is  used  for  location 
calculations  only  if  the  epicenter/station  pair  is  S-P  calibrated. 

In  that  case,  the  arrival-time  error  is  deduced  from  the  P-arrival 
errors,  thus  eliminating  the  need  for  S-wave  travel-time  and  travel- 
time  s.d.  tables.  Since  is  event-dependent,  the  S-wave  file 

name  or  description  should  indicate  the  relevant  event.  For  the 
event,  the  tectonic  media  attenuation  file  is  not  required. 
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1.  Name  and  description  of  the  file, 

2.  Stable  media  attenuation  file. 


*  S-wave  stable  media  attenuation  table  must  indicate 
whether  the  attenuation  table  applies  to  a  shallow  or 
deep  event,  since  the  shallow  attenuation  table  starts 
with  the  Sn  phase. 

3.  Tectonic  media  attenuation  file  (d  $  20”). 

4.  [KEl^Si,  o  -  S  or  T. 

5.  a  -  S  or  T. 

6.  k  -  4. 

7.  ratio  of  P-to-S-wave  velocities  for  S-P  calibrated 
paths . 


pP-Wave  File 

The  pP  wave  is  generated  in  SNAP/D  only  by  a  event. 


1. 

2. 

3. 

4. 

5. 

6. 
7, 


Name  and  description  of  the  file. 
Attenuation  file. 

Travel- time  file. 

Travel-time  s.d.  file. 


iraik-s- 

(ralic,. 


The  foVLauting  three  inpute  are  used  by  the  SUBSET  module  as 
explained  in  Sea,  IV. 


8.  L,  the  total  number  of  step-out  rings. 

9.  L^,  the  mlnlnium  number  of  detecting  rings  required  to  satisfy 
step-out  criteria. 

10.  i  ■  1,  ...,  L,  the  minimum  radius  of  ring  1. 


Ravleigh-Wave  File 

The  Rayleigh  wave  is  not  used  for  location.  The  file  name  or 
description  should  indicate  the  relevant  event  type. 
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IV .  MICROSTRUCTURE 

This  section  sunmarizes  the  calculations  performed  by  each  SNAP/D 
nodule.  We  begin  with  a  review  of  the  different  computational  sequences 
Initiated  by  the  user  during  the  Interactive  session. 

If  the  run  node  Is  detection  only,  then  ATSD  and  LOCATE  are  not  used 
and  the  only  waves  propagated  are  those  Included  In  the  detection  cri¬ 
teria.  If  the  run  mode  Is  location  only,  TRANS  and  BINMB  are  not  used, 
and  all  waves  for  the  selected  event  type  are  propagated  (see  Table  1, 
p.  8),  except  those  excluded  by  the  user  In  the  comnon-wave  Inputs. 

If  the  magnitude  contours  for  a  fixed  probability  threshold  are 
desired,  then  a  binary  search  Is  carried  out  (under  the  direction  of 
BINMB),  In  the  feedback  loop  consisting  of  BINMB,  WAMP,  SPROB,  and 
TRANS  (see  Fig.  1,  p.  11).  If  the  user  requests  probability  contours, 
then  BINMB  Is  not  used. 

Location  calculations  are  performed  either  at  the  user-specified 
magnitude 

magnitude  for  each  epicenter.  Hence,  the  user  should  not  set  m^^  >  0 
If  (1)  the  run  Is  for  location  only;  (2)  probability  contours  are  re¬ 
quested  for  the  detection  output  (l.e.,  output  option  PR  was  chosen); 

(3)  the  P  wave  Is  not  Included  In  the  detection  criteria. 

It  Is  essential  that  the  user  select  commensurate  Input  files. 

The  epicenter/station  calibration  file  requires  cross-referencing  of 
the  station  and  epicenter  Indices  for  the  proper  operation  of  SNAP/D. 
Furthermore,  wave  files  are  generally  relevant  for  a  specific  event 
type  only.  The  user  will  receive  ample  warning  of  those  contingencies 
during  Interactive  prompting. 

Because  of  the  abundance  of  symbols  and  Indices  required  to  de¬ 
scribe  the  mathematical  manipulations  performed  by  SNAP/D,  the  module 
discussions  In  this  section  will  occasionally  share  symbols.  However, 
the  following  Indexing  system  remains  constant: 


,  or.  If  the  user  enters  m^^^  ■>  0,  at  the  detection  threshold 
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i  ■  station  Index, 
j  ■  epicenter  index, 
k  «  wave  index. 

Finally,  we  again  note  that  the  partitioning  of  SNAP/D  into  modules 
is  conceptual.  The  correspondence  of  modules  with  subroutines  in  the 
computer  program  la  discussed  in  Part  II.  Module  descriptions  follow 
below. 


EGRID 

SNAP/D  performs  all  computations  on  an  eplcenter-to-epicenter 

basis.  The  function  of  EGRID  is  to  step  the  program  from  one  epicenter 

to  another  according  to  the  sequence  listed  on  p.  19  of  Sec.  III. 

EGRID  requires  the  inputs  in  the  epicenter  grid  file,  and  the  EGRID 

output  is  (6_. ,  ♦-.)!  the  location  of  epicenter  j  for  j  »  1,  2,  ...,  J. 
EJ  EJ 

GE(»1 

-  t 

Given  the  location  of  the  epicenter  J,  (0gj»  ♦gj)*  station  i, 

(9'^,  GEGM  calculates  the  eplcentral  angular  distance  and 

azimuth  (measxired  positive  clockwise  from  north) ,  from  the  great- 

circle  geometric  relations  [Beyer,  1978]: 


sin  S'gj  cos  sin  A^^^  +  cos  cos  A^^^ 


(1) 


tan(^^  -  ♦gj) 


sin  sin  A^ 


sin  cos  -  cos  sin  A^^j  cos 


After  some  manipulations  of  Eqs.  (1)  and  (2),  we  obtain 

^Glven  a  latitude  6  in  standard  geogrmhic  convention  (l.e., 

6  ■  0  at  the  equator  and  positive  north),  o  ■  v/l  -  6.  The  use  of  the 
polar  coordinate  S'  simplifies  the  expressions  in  the  discussion  of  6E0M. 
The  longitude  ^  is  the  same  in  the  two  conventions. 
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cos 


sin  6*2 j  sin  9*^  cos  -  ♦gj)  +  cos  cos 


and 


tan  a 


ij 


_ Bln  sin  _ 

sin  cos  <|i^  -  cos  Sgj  sin  cos  (♦^  -  ♦gj) 


For  shallow  events  with  £  20  deg»  GECM  also  calculates  , 
the  ratio  of  the  length  of  the  epicenter-to-station  path  in  tectonic 
media  to  the  total  great'^ircle  pathlength  A^^ .  This  calculation  be¬ 
gins  by  determining  the  pathlength  traversed  from  epicenter  to  station 
in  each  regional  grid.  For  every  grid  boundary  longitude  between 

4)„ .  and  ,  the  distance  from  the  epicenter  to  that  longitude  along 
tj  1 

the  propagation  path  is,  from  Eqs.  (1)  and  (2), 


tan  A'  ■ 


tan  (((»’  -  (|)g^)  sin 


tan  (<|i'  -  ♦gj)  cos  ffgj  cos  +  sin  o^j 


(3) 


Likewise,  for  every  grid  boimdary  latitude  S*  <  S'  11  4  ^ 

and  0  >  oL.  if  a..  >  |if/2|,  the  distance  from  the  epicenter  to  that 

EJ  ij 

latitude  along  the  propagation  path  is 


5*  sin  91 


.  cos  V  sin  bp.  cos  a..  ±  cos 
tan  A  - - ^ - 


cos 


9*  cos  91.  -f  sin  91-  cos  a.. 

Ej  Ej 


(4) 


where 


Q  •  sin^  91 


EJ 


2 

cos  a 


iJ 


+  cos 


®EJ  ■ 


cos 


9* 


Since  it  la  possible  to  cross  a  regional  grid  latitude  9*  <  if 
“44  <  K/2|,  or  9*  >  S'  if  a  >  |'n'/2|,  all  real  solutions  of  Eq.  (4) 
represent  admissible  regional  grid  boundary  crossings  unless  A  <  A44 

it 

Now  place  the  {A'}  from  Eq.  (3)  and  the  admissible  {A  }  from 
Eq.  (4)  in  ascending  order  to  form  the  set  A^,  u  1,  ....  U.  Let 
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Aq  =  0  and  ^^2.  ^  ’  Then,  the  nth  path  segment  length  is  6^  ■ 

SQuations  (1)  and  (2)  can  now  be  used  to  find  the  index  of 
the  regional  grid  that  contains  the  nth  path  segment;  let  be  the 
weight  associated  with  this  grid.  All  path  segments  in  the  polar  cap 
regions  (O'  <  AO*^  or  O'  >  it  >  AOj^)  are  ^signed  the  respective  polar 
cap  wel^ts  or  w^.  Then 


Since  0  s  w  i  1  and  7 6  ■  then  0  s  W. .  a:  1.  A  similar  procedure 

u  u  i  j  *  i  j 

is  used  to  determine  the  set  of  Lg  squelch  weights  for  the  epl- 

central  path.  Then 


U+1 

'll 'S  ’u  •  V 

u-1 

since  •  0  or  1,  0  i  i  U  +  1.  If  0,  then  the  eplcentral 

path  is  assumed  to  have  traversed  a  region  that  severely  attenuates  Lg 
waves.  Thus,  if  ^  0,  the  Lg  wave  is  not  propagated  for  that 
epicenter /station  pair  (i.e.,  p^^  -  0)«  If  ■  0,  the  propagation 
of  the  Lg  wave  is  unaffected  by  the  squelch  weights. 

WAMP 

For  a  given  epicenter  J  and  value,  WAMP  calculates  log^^  ^Ijk’ 
the  mean  of  the  logarithmic  amplitude  of  wave  k  observed  at  station  i. 
In  general,^ 

^Although  b|^  and  c^  factors  are  tabulated  by  degrees  in  wave 
attenuation  tables,  the  A, .  In  Eq.  (5)  are  converted  to  kilometers 
In  SNAP/D. 
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loglO 


*  '■f ’«11>  *  ^“*10  ‘ij  ^  'ilk 


(5) 


is  calculated  for  both  stable  (a  -  S)  and  tectonic  (o  ■  T)  media. 
The  event  and  magnitude  factors  [KE]^*^  and  respectively, 

are  used  to  convert  an  m^  value  into  an  m^  value  in  the  magnitude 
units  of  wave  1c.  Thus, 


and,  by  definition,  [KE]^1  -  0  and  -  1.  If  a  single  wave  de¬ 

tection  criterion  is  used,  the  above  expression  is  substituted  into  Eq. 
(5).  The  b^®^(Aj,j)  and  Eq.  (5)  are  the  standard  attenua¬ 
tion  table  entries,  and  the  epicenter/station  calibration  cor¬ 

rection  term  for  wave  k. 

If  wave  k  does  not  require  regional  attenuation,  or  if  >  20*, 
then  logj^Q  is  computed  directly  from  Eq.  (5)  using  the  stable 

media  attenuation  table— the  only  table  available  for  waves  not  re¬ 
gionally  attenuated  and  having  entries  for  A^^j  >  20*  if  the  wave  is 
regionally  attenuated.  Otherwise, 


^°*10  ^Ijk 


(1  ”  W,  )  log..  A...  + 


'10  "ijk  ^  "ij  "''*10  "ijk  • 


where  is  the  regional  path  weight  computed  by  GEOM.  As  discussed 
in  GEOM,  if  the  Lg  squelch  weight  over  the  propagation  path  is  not 
zero,  then,  in  effect,  logj^Q 


SPROB 

SPROB  computes  probability  that  wave  k  propagated  from 
epicenter  j  will  be  observed  at  station  1.  The  following  statistical 
parameters  are  needed  to  compute 
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R^,  the  reliability  of  station  i. 

the  mean  level  of  the  lognormally  distributed  noise 
amplitude  for  wave  Ic  at  station  i. 
the  s.d.  of  the  noise  amplitude  for  wave  k  at 

station  i. 

the  s.d.  of  the  signal  amplitude  for  wave  k  at 

station  i  for  a  given  magnitude  m^,  in  the  magnitude 
units  for  wave  k. 

^bk*  additional  s.d.  of  the  ^og^^  signal  amplitude  for 
wave  k,  given  an  m^  value  (hence,  =  0). 

^Ik*  necessary  for  observation  of  wave  k  at  sta¬ 

tion  1, 

logio  mean  of  the  logarithmic  amplitude  of  wave  k,  gen¬ 

erated  by  epicenter  j  and  observed  at  station  i  (as 
computed  by  WAMP). 

Noting  that  both  signal  and  noise  are  lognormally  distributed, 
and  the  conversion  from  natural  to  base  10  logarithms,  SPROB 

calculates  [Wirth,  1977] 


If  the  user  specifies  a  detection  criterion  Involving  a  single  wave, 
then  detection  contours  are  computed  for  the  magnitude  units  of  that 


wave. 


In  this  case,  =  0  in  Eq,  (6). 
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TRANS 

The  required  notation  for  the  TRANS  and  NPROB  modules  can  be 
Illustrated  by  an  example  of  the  user-entered  detection  criteria  In 
an  Interactive  format: 

(P  *  LG)/2  +  P/4  ,  (7) 


which  literally  means  that  a  network  detection  consists  of  at  least 
a  two-station  observation  of  the  P  cmd  (*)  Lg  waves,  or  (+)  at  least 
a  four-station  observation  of  the  P  wave.  Expression  (7)  Is  said  to 
consist  of  subcriteria  (P  *  LG)/2  and  P/4,  %dilch  specify  wave  combina¬ 
tions  -  P  *  LG  and  2^  ■  P,  respectively. 

The  sequence  of  operations  performed  on  detection  criteria  such 
as  Eq.  (7)  falls  Into  two  categories.  First,  at  the  Initiation  of 
each  run,  the  criteria  are  decomposed  and  reduced  to  a  set  of  canoni¬ 
cal  probabilities.  Second,  at  each  epicenter,  for  a  given  m^  value, 
and  with  the  help  of  the  Individual  wave/station  probabilities  computed 
by  SPROB,  TRANS  calls  on  NPROB  to  perform  the  canonical  calculations. 
Those  results  are  then  used  to  compute  the  network  detection  probabil¬ 
ity,  which  Is  fed  to  BINMB. 

The  sequence  of  TRANS  operations  at  the  start  of  each  SNAP/D  run 
requiring  detection  contours  consists  of  the  following  steps: 


1.  Interactively  Inform  the  user  of  syntax  errors  In  the  user- 
specified  detection  criteria.  The  error  messages  are  covered 
In  Part  II  of  this  report. 

2.  Reorder  the  subcriteria,  the  waves  appearing  In  each  wave 
combination,  and  the  logical  operations,  so  that  the  entire 
expression  Is  in  a  standard  format.  For  example:^ 

^The  examples,  taken  directly  from  TRANS  output,  merely  illus¬ 
trate  the  operation  of  TRANS  and  have  no  operational  significance. 
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a.  PP/3  *  S/1  +  (PP  +  LG)/2  becomes  (LG  +  PP)/2  *  (S/1  *  PP/3). 

b.  (LG  +  P  *  S  +  PP)/4  *  (LG  +  P  *  S)/5  +  (P  *  S  +  PP)/6 
becomes 

(P  *  S  +  LG)/5  *  (P  *  S  +  (LG  +  PP))/4  +  (P  *  S  +  PP)/6. 

3.  Eliminate  logical  subsets  among  reordered  subcriteria  (i.e., 
remove  redundant  subcriteria) ,  and  reduce  wave  combinations 
to  their  simplest  form.  For  example  (from  step  2) : 

a.  (LG  +  PP)/2  *  (S/1  +  PP/3)  becomes  S/1  *  PP/3. 

b.  (P  *  S  +  LG)/5  *  (P  *  S  +  (LG  +  PP))/4  +  (P  *  S  +  PP)/6 
becomes  (P  *  S  +  LG)/5  +  (P  *  S  +  PP)/6. 

4.  Transform  the  reduced  logical  expression  into  an  algebraic 
expression  involving  the  marginal  probabilities  of  independent 
subcriteria  and  the  joint  probabilities  of  dependent  pairs  of 
subcriteria.  This  transformation  thus  eliminates  all  logical 
"o3»s"  among  subcriteria  by  use  of  the  elementary  rule 

+  5^/m]  -  ^[#/n]  +  ^S^/m]  -  *  ^/m]  , 

wtl  t7l2 

where.  If  m  and  m  have  no  waves  In  connion, 

*  i;^/m]  -  . 

For  example  (from  step  3) ; 

a.  f\Sll  *  PP/3]  - 

b.  ?*[(?  *  S  +  LG)/5  +  (P  *  S  +  PP)/6] 

-  ^(P  *  W  LG)/5]  +  *  S  +  PP)/6] 

-  ^(P  *  S  +  LG)/5  *  (P  *  S  +  PP)/6]. 

Interactively  Inform  the  user  of  an  overdependency  among  sub- 
criteria.  Overdependency  with  respect  to  TRANS  exists  if  and 
only  if  a  wave  combination  has  one  or  more  waves  in  common 
with  two  or  more  other  wave  combinations  within  a  probability 
argument  after  step  4.  For  example: 
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Suppose  2r  has  a  wave  in  coiomon  with  both  ST  and  Jr,  and 
step  4  resulted  in  the  expression  *- 

*  ^/n2  *  . 

The  user  would  then  be  informed  of  an  overdependency ,  and 
the  run  would  be  terminated. 

6.  Simplify  the  algebraic  expression  by  adding  and  subtracting 
like  terms.  After  this  step,  the  user  is  interactively  pro¬ 
vided  with  the  resulting  expression. 

7.  Obtain  and  algebraically  reduce  the  expression  for  the  prob¬ 
ability  of  observing  the  wave  combination  at  a  given  station 
for  each  Independent  wave  combination.  Likewise,  for  each 
dependent  pair  of  wave  combinations,  obtain  and  reduce  ex¬ 
pressions  for  the  marginal  probabilities  of  observing  each 
wave  combination,  and  the  probability  of  observing  both  com¬ 
binations  jointly.  For  example  (from  step  5,  example  b) : 

a.  The  expression  is  of  the  form 

*  :0^/6]  . 

The  probability  of  observing  at  station  1  is 

^i  ”  ^11^14  ^  ^13  ■  Pil**  13^14  • 

Hie  probability  of  observing  m  at  station  i  is 

t 

\  “  PilPi4  Pi5  "  ^11^14^15  * 

X  2 

The  probability  of  observing  2  *  2  at  station  1  is, 

after  reduction, 

*^1  “  PllPl4  Pl3*’i5  '■  PilPl3*’i4Pi5  * 

The  p^j^  in  the  above  expressions  are  the  probabilities 
of  observing  wave  k  at  station  i  for  a  given  epicenter 
(we  have  suppressed  the  j)  and  m^  value. 

This  concludes  the  operations  performed  by  TRANS  at  the  start  of  a  run 
requiring  the  computation  of  detection  contours. 
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The  next  sequence  of  operations  performed  by  TRANS  Is  done  for 
each  m^  value  at  every  epicenter.  The  sequence  begins  with  the  re¬ 
ception  of  i  •  1|  ....  N,  k  ■  1,  . . . ,  K,  from  SPROB.  It  con¬ 

cludes  with  the  computation  of  the  network  detection  probability  P^. 
As  an  intermediate  step,  TRANS  calls  on  NPROB  to  compute  the  marginal 
probability  of  Individual  subcriteria  and  the  joint  probability  of 
dependent  pairs  of  subcriteria.  The  sequence  proceeds  as  follows: 


8.  Compute  the  numerical  value  of  {p^},  1  >  1,  ...,  N,  asso¬ 
ciated  with  each  individual  term  of  the  form  and 

{p^,  q^,  1  ~  1,  ...,  N,  associated  with  each  joint  term 

of  the  form  *^/m],  using  from  SPROB  and  the 

probability  expressions  from  step  7. 

9.  Initiate  a  one-dimensional  (ID)  calculation  for  each  indi¬ 
vidual  term  by  supplying  NPROB  with  the  {p^^}  associated  with 
that  term.  NPROB  returns 

probability  that  exactly  n'  out  of  N  stations  observe  wave 
combination  3.  Clearly, 

N 

n'*n 

which  is  computed  by  TRANS  for  each  individual  term  in  the 
expression  resulting  from  step  6. 

10.  Initiate  a  two-dimensional  (ZD)  calculation  for  each  joint 

term  by  supplying  NPROB  with  the  {p^,  q^,  c^}  associated  with 
that  term.  NPROB  returns  the  proba- 

blllty  that  exactly  n'  out  of  N,  and  m*  out  of  N,  stations 
observed  both  ^  and  respectively.  Clearly, 

n'*n  m'"m 

which  is  computed  by  TRANS  for  each  joint  term  in  the  expres¬ 
sion  resulting  from  step  6. 
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11.  Calculate  from  the  results  of  steps  6,  9,  and  10. 

NPROB 

Each  call  to  NPROB  produces  one  of  two  types  of  canonical  prob¬ 
ability  calculations  for  TRANS.  The  calculations,  many  of  which  may 
be  required  for  each  m^  value  at  a  given  epicenter,  include  general¬ 
izations  of  ID  and  2D  discrete  convolutions  [Feller,  1968].  The  ID 
calculation  is  essentially  the  same  procedure  used  in  NETWORTH  [Wlrth, 
1977];  the  2D  calculation  was  developed  to  accommodate  multiwave  net¬ 
work  detection  criteria  that  use  combinations  of  dependent  wave 
arrivals  at  individual  stations,  as  explained  in  the  discussion  of 
TRANS.  The  2D  convolution,  performed  %rlthout  approximation,  is  the 
most  computationally  efficient  procedure  among  the  several  tested  for 
the  2D  calculation. 

ID  Calculation 

Provided  by  TRANS  with  the  probabilities  of  observing  a  given 
wave  combination^  at  each  station,  {p^},  1*1,  ...,  N,  NPROB  calcu¬ 
lates  ^he  probability  of  observing  3  at  exactly  n  of  the 

total  N  stations  in  the  network.  Since  here  we  always  refer  to  wave 
combination  3,  we  generally  denote  the  probability  of  observing  3  at 
exactly  n  of  M  s  N  stations  with  indices  Yj^t  •••»  Yjj  by  P(n|Yj^,  ...,  Yj^ 
so  that  *"» 

The  general  ID  convolution  equation  is  given  by 

“2 

•••*  Yj^)  “  ...,  )F(n  —  u|82^,  .••»  )  , 

U*Uj^  ^  ^ 

where  M  ■  +  M2  , 

-  max(0,  n  -  M2)  , 

U2  ■  mln(n,  Mj^)  . 
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Note  that  the  division  of  the  M  stations  into  the  two  sets  of  and 
M2  stations  Is  arbitrary  If,  and  only  If,  the  two  sets  are  mutually 
exclusive  and  exhaustive  with  respect  to  the  M-statlon  set.  We 
schematically  represent  the  combining  of  station  sets  In  the  above 
expression  by 


(Mj^)  (M2) 


(M)  . 


In  order  to  show  how  NPROB  calculates  P(n|l,  ....  N),  consider 
an  N  ■  8  station  network.  At  the  first  stage,  calculate 


P(n|i,  1  +  1)  P(n|i)P(n  -  u|i  +  1)  ,  n  -  0,  1,  2  , 

U«Uj^ 

for  1  -  1,  3,  5,  7.  Note  that  P(l|i)  -  p^  and  P(0|i)  -  1  -  p^, 
1*1,  ...,  N.  At  the  second  stage,  calculate 


P(n|l,  ...,  A)  P(n|l,  2)P(n  -  ul3.  A)  ,  n  -  0,  1,  ...,  A 

U«Uj^ 

and  similarly  obtain  P(n|5,  ....  8).  Finally,  at  the  third  stage, 
calculate 


PCnjl,  ....  8) 


P(n|l»  A)P(n  -  ul5, 

U«Uj^ 


. . . ,  8)  , 
n  •  0,  1,  ...,  8 


This  example  is  simplified  by  the  fact  that  N  ■  8  is  a  power  of 
two.  Thus,  the  combining  rule  is  represented  by 
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Stage 


1 

2 

3 


Combinations  of  Stations 

(1)  (1)  (1)  (1)  (1)  (1)  (1) _ a) 

(2)  (2)  (2)  (2) 

(4)  (4) 

(8) 


However,  the  general  ID  convolution  also  allows  for  the  combining  of 
different  sized  sets  of  stations.  Thus,  NPROB  would  handle  an  N  ■  23 
station  network  as  follows: 


Stage 

1 

2 

3 

4 

5 


Combinations  of  Stations 
(1)  (1)  (1)  (1)  •••  (1)  (1)  (1) 
(2)  (2)  •••  t2)  (1) 

(4)  (4)  (4)  (4)  (4)  (3) 


(8)  (8)  (7) 


(23) 


2D  Calculation 

TRAHS  initiates  a  2D  calculation  by  providing  NPROB  with  the  proba¬ 
bilities  {p^,  q^,  c^},  i  <■  1,  ...,  N,  as  defined  in  steps  7  and  10  of 
the  TRANS  subsection.  NPROB  then  computes 

Joint  probability  of  observing  wave  combinations  and  ^  at  exactly 
n  and  m  stations,  respectively,  of  the  total  of  N  stations.  With 
and  ^  understood,  we  generally  denote  the  probability  of  observing 
at  n,  and  at  m,  of  the  M  $  N  stations  with  indices  •••• 
as  P(n,  bIYj^,  ....  Yjj).  Thus, 

^^(n:N)  *^m:N)^  “  * 
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The  2D  calculation  begins  with  the  computation  of  P(n,  m|l), 

1  ~  1,  N.  We  will  need  the  following  elementary  set  relations 

for  events  and  G2: 


Gj^G2  ”  ^  ^  Morgan's  Law) 


and 


G,  -  G,G2  • 

Letting  and  ^q.^)  “  ^o*'  station  1,  we  obtain 


P(n,  m|l)  - 


1  -  Pi  -  . 

. 


9.  -  c. 


■1  ’ 


n  ■  m  ■  0  , 
n  ■  1,  m  ■  0  , 
n  «  0,  m  “  1  , 
n  •  m  ■  1  , 


for  1  ■  1,  . . . ,  N. 

The  general  2D  convolution  equation  for  combining  station  sets  of 
size  and  M  ■  +  M2»  Is 

P(n,  mjy^,  •••!  Yjj) 


where  u^^  ■  max(0,  n  -  M2), 
U2  -  mln(n,  Mj^), 

-  max(0,  m  -  M2), 
V2  ■  mln(m,  Mj^). 


“2  "^2 


^i“i’  °Mj^^ 


U»Uj^  V^J^ 


X  P(n  -  u,  m  -  v|6j^,  ...,  Bjj  )  , 

n,  m  ■  0,  1,  . . . ,  M  , 
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Using  this  expression  and  P(n,  n|i)  from  above,  NPROB  combines 
station  sets  in  the  same  manner  as  the  ID  calculation.  For  example, 
at  stage  1, 


P(n,  m|i,  i  +  1) 


“2  ^2 

P(u,  v|i)P(B  -  u,  m  -  vjl  +  1)  t 


n“0,  1,  2  , 


for  i  1,  3,  5 .  Then,  for  stage  2, 

“2  ^^2 

P(n,  m|l,  2,  3,  4)  -  ^  ^  ^  ^  P(u,  vjl,  2)P(n  -  u,  m  -  vl3,  4)  , 

X1*0^  Xy  esej  A  f 

etc.,  until  P(n,  m|l,  ...,  N),  n,  m  ■  0,  1,  ...,  N,  is  obtained.  The 
general  combining  scheme  for  arbitrary  N  is  the  same  as  that  illus¬ 
trated  for  the  ID  calculation. 

BINMB 

BINMB  directs  the  binary  search  for  the  magnitude  value  (gen¬ 
erally  in  units)  that  results  in  P^,  the  user-specified  threshold 
probability  for  network  detection.  Initially,  BINMB  selects  m^  ■  ^mln* 
the  user-specified  minimum  magnitude  for  the  search.  WAMP,  SPROB,  and 
TRANS  then  perform  the  calculations  necessary  to  obtain  Fj^(m^),  the 
network  detection  probability  for  magnitude  m^.  If  then 

BINMB  selects  m2  ■  and  receives  TRANS.  If  ^  ^T* 
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If  PQ(m3)  <  P^,  then 


m  +  m 
"4 - — 


Otherwise, 


m. 


+«3 


etc.,  until 


|P(«^)  -  P^l  ^  0.01 

for  some  v.  The  threshold  magnitude  is  then  m^,  and  the  search  is 
terminated. 

ATSD 

ATSD  computes  the  arrival-time  s.d.  for  all  waves  that  provide 
the  arrival-time  data  that  LOCATE  uses  for  computing  the  precision  of 
the  hypocenter  estimate.  Those  waves  Include  the  P  (k  ■  1),  Pg  (k  ■  2), 
S  (k  ■  4),  and  pP  (k  ■  5)  waves.  The  Lg  wave  (k  ■  3)  is  also  used  for 
location.  However,  Lg  observations  provide  back  azimuth  data  rather 
than  arrival-time  data.  The  Rayleigh  wave  (k  >  6)  is  not  used  for 
location. 

The  general  expression  for  computing  the  arrival- time  s.d.  for 
wave  k,  propagated  fr'oin  epicenter  j  and  observed  at  station  i,  is 

"uk  ■  K + 4  * 

where  is  the  reading-error  s.d.,  is  the  station-dependent 
arrival-time  s.d.,  and  the  travel-time  s.d.  table  entry  for 

wave  k.  The  the  Pg  and  pP  waves  are  computed  from  Eq.  (8) 

without  modification.  However,  the  P  and  S  waves  require  special 
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treatment  because  of  the  possibility  of  epicenter/station  master-event 
or  S-P  calibration. 

Recall  that  is  the  set  of  indices  for  stations  that  are  master- 
event  calibrated  with  respect  to  epicenter  J .  The  effect  of  calibrating 
the  eplcenter-to-statlon  path  is  to  reduce  scatter  in  the  P-wave  arrival 
time  [Clervo,  1978].  However,  the  reading  error  is  always  present. 

Thus,  for  the  P  wave,  ATSD  computes  if  i  €  ;  otherwise, 

is  computed  from  Eq.  (8). 

The  situation  with  respect  to  the  S  wave  is  more  complex.  The 
S  wave  is  used  for  location  only  if  station  i  is  S-P  calibrated  with 
respect  to  epicenter  J  (l.e.,  if  1  €  Thus,  if  i  £  ,  then  LOCATE 

receives  =  0,  effectively  eliminating  the  S  wave  from  location 

calculations.  If  1  €  ({>j,  then  LOCATE  considers  the  error  in  the  datum 
to  be  composed  of  a  linear  combination  of  the  P-arrival  time  (t^)  and 
the  S-arrival  time  (tg)  given  by 


This  particular  combination  of  tp  and  tg  is  uncorrelated  with  respect 
to  tp,  which  simplifies  the  matrix  algebra  in  LOCATE.  The  s.d.  of  T 
is 


Thus,  if  i  €  then  "  '^Ol  discussed  above,  and  -  Oq^, 

the  S-wave  reading  error.  In  the  SHAP/D  computer  program,  LOCATE  is 
fed  and  Vp/vg  at  the  beginning  of  the  run.  Then,  is 

computed  within  LOCATE  for  each  (i,  J)  from  Eq.  (9)  using  the  value 
of  supplied  by  ATSD. 
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SUBSET 

The  function  of  SUBSET  is  to  perform  all  computations  involving 
subsets  of  the  entire  N-station  network.  In  SMAP/D,  these  Include 
calibrations  concerning  the  step-out  phenomena  for  the  reflected  phase 
pP.  In  SIIAP/I  (planned  for  future  development),  SUBSET  would  also 
handle  calculations  for  the  "first-motion"  discriminant  [Ciervo,  1978]. 

The  step-out  requirement  specifies  that  the  pP  wave  must  be  ob¬ 
served  over  a  range  of  epicentral  distances  in  order  to  be  correctly 
identified  [Ciervo,  1978],  The  identification  process  is  siiiq>ly  to 
observe  that  the  pP  arrival  lags  farther  behind  the  P  arrival  as  epi¬ 
central  distance  increases  (hence  the  term  "step-out").  As  shown  in 
Fig.  2,  SNAP/D  models  the  step-out  process  by  calculating  the  proba¬ 
bility  that  the  pP  wave  will  be  observed  by  at  least  one  station  in 
of  a  total  of  L  "rings"  surrounding  the  epicenter.  L^,  L,  and 
i  •  1,  L,  the  minimum  radius  of  ring  £.,  are  all  entered  by  the 

user  into  the  pP-wave  file.  Stations  for  which  <  A^  are  assumed 
not  to  observe  the  pP  wave. 

Station  i  is  said  to  be  in  ring  £  if  A^^  s;  A^j  <  where  A^^j^ 

is,  by  definition,  the  distance  corresponding  to  the  last  pP-wave 
attenuation  table  entry.  Furthermore,  ring  £  is  said  to  observe  the 
pP  wave  if  at  least  one  station  in  ring  £  observes  it.  Now  define  the 
events 


pP  propagated  from  epicenter  j  is  observed  by  station  1  , 


and 


Lg  of  L  rings  with  respect  to  epicenter  j  observe  the  pP  wave  . 


As  computed  by  SPROB, 
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The  function  of  SUBSET  in  SNAP/D  Is  to  calculate 


(11) 


In  order  to  calculate  we  will  also  need  to  define  4^ 

as  the  set  of  indices  for  all  stations  in  ring  Then  the  probability 
that  ring  I  observes  the  pP  wave  is 


1-1, 


L  . 


The  conditioning  in  the  expression  ensures  that  the  ring 

containing  station  i--call  it  ring  1’ —observed  the  pP  wave.  Thus, 

is  the  probability  of  L-  -  1  "successes"  out  of  L  -  1  "trials," 
where  the  L  -  1  trials  have  the  success  probabilities  ,  . . . ,  , 

excluding  Q^, .  This  probability  is  Identical  to  that  computed  by  the 
ID  calculation  in  NFROB.  Thus,  SUBSET  calls  NPROB  with  the  set  {Q£)> 
1-1,  . . . ,  L,  excluding  1* ,  and  receives  back 


as  the  probability  of  exactly  Lg  -  1  successes  in  L  -  1  trials.  Then 
is  calculated  from  Eq.  (11).  Por  a  given  epicenter  j,  the 
i  -  1,  ...,  N,  are  then  used  instead  of  the  1  “  1.  N,  for 

the  computations  in  TKANS  and  LOCATE. 


LOCATE 

For  epicenter  j ,  LOCATE  calculates  A^ ,  the  area  of  the  eplcentral 
estimate  uncertainty  ellipse;  and  d^,  the  length  of  the  depth'-estimate 
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confidence  interval  at  P  ,  the  user-specified  confidence  probability. 

c 

For  each  set  of  epicenter  indices,  T  ,  m  ■  1,  ...»  M,  LOCATE  also 

n 

calculates  the  minimum,  maximum,  and  average  Importance  values  for 
each  station/wave  pair.  The  inputs  required  by  LOCATE  for  each  epi¬ 
center  are  ^°ijk^*  ^‘^$1^’  ^^m^’  ^®lj^’  ^P’  ^c* 

and  the  travel-time  tables  for  the  relevant  waves. 

LOCATE,  a  complex  computer  model  by  itself,  was  developed  by 
Systems,  Science  and  Software,  which  supplied  the  general  description 
provided  in  the  Appendix, 

PRINT 

The  fimctlon  of  PRINT  is  to  format  the  output  of  SNAP/D.  The 
types  of  output  are  determined  by  the  user  during  the  interactive 
session  at  the  initiation  of  each  run.  Since  much«of  the  useful  out¬ 
put  for  each  mn  is  the  input  material  itself,  the  user  may  wish  to 
have  a  printout  of  the  entire  CRT  interactive  session.  The  user  may 
also  elect  to  have  all  table  (e.g.,  travel  time  and  attenuation) 
entries— material  that  would  not  be  seen  if  input  files  were  called 
without  review  or  update — printed  out  for  any  selected  wave.  In  addi¬ 
tion,  summary  information  or  detailed  results  may  be  printed  out  for 
each  epicenter. 

SNAP/D  output  consists  of  three  sections.  The  first  section 
identifies  the  run  mode  (detection,  location,  or  both) ;  provides  the 
network  probability  expression;  and,  depending  on  the  mode,  gives  the 
threshold  magnitude,  probability,  or  percent  confidence  level  for  the 
ellipses.  These  are  followed  by  the  user's  file  descriptions  of  every 
file  used. 

The  second  section  consists  of  station  and  wave  tables.  The  sta¬ 
tion  tables  are  printed  one  table  per  wave.  Each  contains  the  station 
index,  name,  and  location;  and  the  station-noise  amplitude  and  s.d., 
followed  by  the  signal-to-nolse  ratio,  arrival-time  s.d.,  station  re¬ 
liability,  and  signal  s.d. 
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Next,  the  set  of  wave  tables,  if  requested,  Is  printed  by  wave 
type.  Three  types  of  table  are  included:  (1)  stable  and  tectonic 
(if  applicable)  attenuation  tables,  which  list  angular  distance  in 
degrees,  b>f actors,  and  c-f actors;  (2)  the  travel-time  table,  lAich 
contains  distance  in  degrees  followed  by  the  number  of  seconds  needed 
to  traverse  that  distance;  (3)  the  travel-time  s.d.  table,  which  lists 
angular  distances  in  degrees,  followed  by  travel-time  standard  devia¬ 
tions  in  seconds. 

The  third  section  (optional)  consists  of  the  SNAP/D  calculations 
for  each  epicenter.  The  first  line  of  the  detailed  epicenter  output 
lists  epicenter  number,  location,  and  threshold  magnitude  for  network 
detection.  The  second  line  lists  station-specific  data  for  each  epi¬ 
center.  Those  consist  of  station  index  and  name,  angular  distance, 
and  azimuth.  Also  Included  for  each  wave  used  is  the  s.d.  of  arrival 
time  (zero  if  the  run  is  for  detection  only),  and  probability  of  ob¬ 
servation  (P^jjj)  for  the  given  station,  epicenter,  and  wave  combination. 

Regardless  of  whether  a  detailed  epicenter  output  is  elected,  a 
network-detection  summary  table  will  be  printed  if  the  run  is  for  de¬ 
tection  only,  detection  and  location  with  m^^  #  0,  or  detection  and 
location  with  m^^^  «  0  and  probability  contours  requested.  The  heading 
of  this  table  lists  the  threshold  magnitude  or  probability,  and  the 
binary-search  magnitude  range  for  a  threshold  magnitude.  Table  entries 
consist  of  epicenter  number,  location,  and  threshold  magnitude  or 
probability . 

For  a  location  or  detection/ location  run,  a  network- location  sum¬ 
mary  table  is  presented  next.  The  heading  of  such  a  table  shows  the 
binary-search  magnitude  range  if  m^^  >  0.  The  columns  of  the  table 
consist  of  epicenter  index,  location,  detection  magnitude  (equal  to 

or  threshold  magnitudes  if  m^^  >  0),  the  lengths  of  the  semi- 
major  and  semifflinor  ellipse  axes,  strike,  ellipse  area,  and  depth- 
interval  length. 

If  importances  are  calculated,  those  tables  are  presented  last. 

They  are  headed  by  importance-set  index  and  consist  of  station  index 
and  the  minimum,  average,  and  maximum  importance  values  for  each 
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station/wave  used  for  location.  The  first  row  of  information  for  each 
station  is  the  location-ellipse  importance  value,  and  the  second  is 
the  depth- importance  value.  If  the  minimum  is  100  and  the  average  and 
maximum  are  both  zero,  the  station  was  not  used  for  location  and, 
hence,  not  used  in  the  importance  section. 


1 

i 


-48- 


1 


AC90C114 


PART  II 
CODE  DOCUMENTATION 
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V.  STRUCTURE  OF  THE  PROGRAM 


In  Part  I,  the  analytical  structure  of  the  SNAP/D  computer  model 
was  presented.  Part  II  details  the  pvogvam  modules,  which  are  simply 
groups  of  subroutines.  The  conceptual  modules  in  Part  I  were  written 
using  all  capital  letters.  In  Part  II,  conceptual  modules  are  indi¬ 
cated  by  capital  italics,  whereas  subroutine  names  are  in  capital 
letters. 

INPUT  STRUCTURE 

Entry  and  exit  from  the  SNAP/D  software  system  occurs  in  the 
"root"  program  module  SNAPD.  That  routine  provides  Initialization  of 
the  computational  environment  and  selection  from  the  main  menu. 

Five  major  routines  called  from  SNAPD  are;  (1)  CALC,  which  in¬ 
terfaces  with  the  set  of  calculative  routines  beginning  with  MAINX; 

(2)  CTLDF,  which  provides  menu  selection  support  for  control  inputs; 

(3)  DETDF,  i^lch  provides  menu  selection  and  data  base  I/O  support  for 
detection  inputs;  (4)  SIGDF,  which  provides  menu  selection  and  data 
base  I/O  support  for  all  wave-related  Inputs;  (5)  TOPODF,  which  pro¬ 
vides  menu  selection  and  data  base  I/O  support  for  all  topographic 
Inputs . 

In  support  of  data  base  I/O  are  routines  ISEL,  lOPEN,  and  ICLOS. 
CTLDF,  DETDF,  SIGDF,  and  TOPODF  call  on  "worker"  routines  to  provide 
interactive  dialogue,  data  base  initialization,  and  data  base  orga¬ 
nization,  so  that  the  data  are  acceptable  to  the  calculative  routines. 
CTLDF  calls  on  CTLPDU,  and  DETDF  calls  on  DETDU.  SIGDF  calls  on  WAVDU 
for  nonspecific  wave  input,  and  on  PWADU,  PGWDU,  LGWDU,  SWDU,  PPWDU, 
and  RAWDU  for  wave  Inputs  specific  to  P,  Pg,  Lg,  S,  pP,  and  Rayleigh 
waves,  respectively.  TOPODF  calls  on  STADU,  EPGDU,  ESCDU,  REGDU,  and 
IMMDU  for  information  specific  to  station  data,  epicenter  grid  data, 
epicenter /station  calibration  data,  regional  grid  data,  and  importance- 
matrix  data,  respectively. 
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MATDU  supports  Interactive  dialogue  pertaining  to  stable  and 
tectonic  media  Information,  and  TTIDU  and  TTSDU  support  Interactive 
dialogue  for  travel-time  data  and  travel-time  standard  deviations  for 
routine,  consists  of  INPUTC  calling  for  a  one-character  Input  (used 
and  INPDLO  provide  a  common  facility  for  entering  latitude  and  longi¬ 
tude  Information  for  the  topographic  work  routines. 

Terminal  Input  support,  which  requires  Interfacing  with  the  LOGGER 
routines,  consists  of  INFUTC  calling  for  a  one-character  Input  (used 
mostly  for  terminal  and  terminal  operator  synchronization) ;  INPUTD 
calling  for  a  one-character  Input  while  prompting  a  default  character; 
INPUTF  calling  for  a  file-name  Input;  INPUTM  calling  for  a  numeric¬ 
string  Input;  and  INPUTS  calling  for  a  character-string  Input.  Routine 
INPUTN  calls  on  CHECKN  for  the  numeric-string  validation.  A  SNAP/D 
Interactive  section  calls  on  STRF  for  a  variety  of  character-string 
manipulation  functions. 

INTERNAL  CODE  STRUCTURE 

Table  2  lists  the  conceptual  modules  and  the  particular  subroutines 
used  In  Implementing  them.  TRANS  Is  performed  by  many  subroutines, 
most  with  a  specific  algebraic  function.  Subroutine  TRANS  is  the  main 
driver  for  breaking  the  detection  criterion  Into  its  subcriteria.  TRANS 
converts  the  Input  detection  criteria  Into  a  form  that  can  be  evaluated. 
It  calls  RULES,  which  directs  the  relevant  algebraic  manipulative  rou¬ 
tines,  Including  SYNTAX,  which  checks  for  syntax  errors  In  the  sub- 
criteria  string;  CHECK,  which  determines  if  the  string  Is  valid;  and 
WRITE,  which  writes  the  string  to  the  terminal.  POUT,  also  called  by 
RULES,  writes  out  Intermediate  forms  of  the  criterion  string. 

The  functions  of  EGRID,  BINMB,  and  WAMP  are  contained  entirely  In 
MAINX.  SPROB  Is  performed  by  MAINX  with  the  probability  function  being 
calculated  In  subroutine  PROBF.  GEOM  Is  performed  by  BNDRY  and  IFIND. 
BNDRY  calculates  the  angular  distances  and  azimuths  between  epicenter 
and  stations.  It  also  finds  the  distance  traveled  in  each  regional 
grid  and  calls  IFIND.  IFIND  determines  the  corresponding  regional  grid 
Index  for  each  section  of  the  epicenter/station  path. 
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Table  2 

CORRESPONDENCE  OF  MODULES  WITH  SUBROUTINES 


Module 

Subroutine 

Name 

TRANS 

TRANS 

RULES 

ADD 

ALTOPO 

ASSOC 

CHECK 

CCRMUT 

DISASS 

DISTRB 

NEGOUT 

POTOAL 

POUT 

PROBAB 

REPLOP 

SPROB 

SUBS 

SVBSET 

WRITE 

SYNTAE 

EGRID 

MAINE 

BINMB 

MAINE 

WAMP 

MAINE 

SPROB 

MAINE 

PROBF 

GEOM 

BNDRY 

IFIMD 

NPROB 

NUMBER 

QPROB 

NETPRB 

PRBCON 

PROBK 

ATSB 

ATSD 

SUBSET 

SUBSET 

PRINT 

LOCATE^ 

MAINE 

PREPRT 

0 

Described  by  Systems,  Science  and 
Software. 


The  rest  of  the  subroutines  are  algebraic.  ALTOPO  converts  from 
algebraic  form  to  reverse  polish,  whereas  POTOAL  converts  from  reverse 
polish  to  algebraic.  ASSOC  performs  the  associative  law;  DISASS  does 
the  opposite.  COMMUT  conmutes  elements  of  relevant  operators  such  as 
addition..  DISTRB  changes  an  event  description  of  the  form  (A  U  B)  0  C 
to  the  form  (A  n  C)  (j  (B  n  C).  ADD  performs  the  addition  operator, 
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whereas  NEGOUT  replaces  negation  with  subtraction.  PROBAB  passes 
the  probability  operator  through  the  criterion,  whereas  SPROB  passes 
the  probability  operator  through  Individual  wave  combinations  to  ob¬ 
tain  the  probabilities  needed  by  the  NPROB  module.  REPLOP  replaces 
algebraic  operators  with  logical  operators.  Finally,  SUBS  finds  and 
eliminates  subsets  within  a  wave  combination,  and  SVBSET  finds  and 
eliminates  subsets  within  subcriteria. 

Module  NPROB  compromises  five  subroutines.  NUMBER  takes  the  out¬ 
put  detection  criterion  from  TRANS^  accepts  probabilities  from  SPROB, 
and  performs  the  network  detection  calculations  via  the  function  QPROB. 
QPROB  then  takes  the  modified  SPROB  output  and  calls  appropriate 
routines  based  on  the  number  of  dependencies.  If  there  is  one  de¬ 
pendency,  QPROB  calls  NETPRB  fWlrth,  1977],  which  finds  the  probability 
of  detection  using  the  function  PROBK.  PROBK  determines  the  proba¬ 
bility  of  finding  an  exact  nuodser  of  events  from  two  groups.  If  there 
are  two  dependencies,  QPROB  takes  the  modified  SPROB  output  and  readies 
it  for  the  2D  convolution.  That  convolution  is  performed  by  PRBCON, 
which  returns  a  final  value  to  QPROB,  yielding  the  detection  probability 
for  the  given  TRANS  output. 

The  ATSD  module  is  wholly  performed  in  subroutine  ATSD.  Module 
SUBSET  Is  performed  by  subroutine  SUBSET,  with  the  conditional  proba¬ 
bility  calculated  by  NETPRB  and  PROBK. 

The  PRINT  module  Is  performed  by  both  MAINE  and  PREPRT.  PREPRT 
writes  out  the  run  mode,  station  tables,  and  desired  wave  tables. 

MAINE  prints  the  summary  network  results,  and,  if  desired,  the  de¬ 
tailed  epicenter  Information. 

The  LOCATE  module  consists  of  routines  supplied  by  Systems,  Science 
and  Software,  which  are  described  In  the  Appendix.  LOCATE  Indicates 
certain  errors,  such  as  a  value  being  outside  table  bounds.  The  module 
Is  fed  Information  from  LOCLOD, 'which  uses  the  output  from  ATSD,  SPROB, 
relevant  Input  parameters,  and  In  the  case  of  a  deep  event,  SUBSET. 
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VI.  INTERACTIVE  AND  BATCH  INPUTS 


GENERAL  COMMENTS 

The  SNAP/D  software  system  was  designed  for  interactive  usage, 
with  batch  processing  as  a  special  case.  Batch  input  consists  of  a 
stream  of  terminal  comnands,  stored  as  a  sequential  file  and  directed 
at  the  SNAF/0  terminal  input  logical  unit  (LUN) .  Further  discussion 
of  these  input  follows  In  the  "Terminal  I/O"  subsection. 

Batch  usage  differs  from  interactive  usage  In  both  LUN  assign¬ 
ment  and  handling  of  pause  ccnanands.  The  following  lists  show  how 
the  LUN  assignments  differ: 


Interactive  LUNs 

Batch  LUNs 

1.  Log  file  (output) 

1. 

Log  print  file 

2.  Data  base  Input/output 

2. 

Data  base  input/output 

3.  Printer  file  (output) 

3. 

Printer  file  (output) 

4 .  Terminal  input 

4. 

Pseudo- terminal  Input 

5 .  Terminal  output 

5. 

Log  print  file 

Pause  commands  are  Ignored  under 

batch 

usage.  The  batch  execution 

halts  when  an  end-of-file  condition  is 

input . 

A  numeric  string  Is  defined  as: 

detected  in  the  pseudo-terminal 

slgn:dlgltl(s) :declmal-polnt:dlgitl(s) :E  or  D:sign:dlglt2(s)  , 

where  the  count  of  digltl(s),  the  precision,  may  not  exceed  15  digits; 
and  the  count  of  dlgit2(s),  the  exponent,  may  not  exceed  2  digits. 

Some  correct  examples  of  numeric  strings  with  similar  values  are 

•  123456789,  +123456789,  1.23456789E+8 

•  -6E5,  -600000,  -6+5 

•  6D-5,  0.000006,  6-5 
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Incorrect  numeric  strings  are  as  follows: 

•  1234567890123456  (precision  Is  too  large) . 

•  -H>5  (multiple  signs). 

•  6D555  (exponent  Is  too  large) . 

•  -£+2  (no  fraction  part). 

•  5E  (partial  exponent  part). 

•  6x5  (unrecognized  character  x) . 

A  numeric  string  can  be  entered  "format-free"  as  long  as  the  string 
Is  contained  within  the  first  20  columns.  The  numeric  string  Is 
terminated  by  the  twenty-first  column,  or  the  first  space  after  the 
start  of  the  numeric  string. 

A  file  name  has  two  parts:  (1)  a  name  consisting  of  eight  or 
fewer  alphanumeric  characters  followed  by  a  period;  (2)  an  extension 
of  three  or  fewer  alphanumeric  characters.  Both  the  name  and  exten¬ 
sion  must  begin  with  an  alphabetic  character.  Table  3  presents  a  list 
of  extension  names  suggested  for  standardized  use.  For  example.  In 
a  particular  seismic  study,  a  file  associated  with  that  study  might 
have  a  file  name  of  STUDY025.X1.  As  can  be  seen  from  Table  3,  this 
would  Indicate  F-wave  data  that  comes  from  an  explosion,  within  the 
aggregate  of  Information  associated  with  the  particular  study  STUDY025. 

TEBMIMAL  I/O  , 

All  terminal  Inputs  are  preceded  by  informational  messages  and 
prompts  that  provide  entry  cues,  such  as  a  name  (SI^IA  004),  a  range  ^ 
(0-1),  or  a  default  value.  The  type  of  reply  Is  Indicated  by  the 
number  of  trailing  periods  on  the  prompt: 

Number  of  Periods  Requested  Terminal  Input 

2  One-character  Input 

3  Numeric-string  Input 

4  Alphanumeric-string  Input 

Default  values  are  selected  by  typing  the  carriage  return  (CR)  key. 


-56- 


AC90C114 


Table  3 

FILE-NAME  EXTENSIONS 


Extension^ 

File  Contents 

Data  Files 

DET 

Detection  Input  data 

STA 

Station  data 

ESC 

Epicenter/station  calibration  data 

EPG 

Epicenter  grid  data 

REG 

Regional  grid  data 

IMM 

Importance-matrix  data 

WAV 

Case  study  wave  file  list 

Seismic  Wave  Data  Files 

Xk 

Wave  data  used  for  an  explosion 

SQk 

Wave  data  used  for  a  shallow  earthquake 

DQk 

Wave  data  used  for  a  deep  earthquate 

ASk 

Stable  media  attenuation  data 

ATk 

Tectonic  media  attenuation  data 

TSk 

Travel-time  data  for  a  shallow  event 

TDk 

Travel-time  data  for  a  deep  event 

SSk 

Travel-time  standard  deviation  for  a 

shallow  event 

SLk 

Travel-time  standard  deviation  for  a 

deep  event 

The  variable  k  denotes  a  wave  type  as  follows: 


k  Wave 

1  P 

2  Pg 

3  . 

4  S 

5  pP 

6  Rayleigh 
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For  a  two-  or  three-period  prompt,  typing  the  "x"  key  will  either 
terminate  a  "logical"  group  of  questions,  or  exit  a  routine  to  a  menu. 
Note  that  typing  the  "x"  key  after  a  four-period  prompt  will  have  no 
effect  other  than  entering  the  character  string  "x,"  and  that  typing 
a  finite  sequence  of  "x"  causes  SHAPD  to  exit  to  the  monitor  level 
of  cononands. 

Three  other  special  keys  may  be  used  with  a  two-  or  three-period 
prompt: 

1.  ;  comments,  which  provide  documentation  (of  terminal  Input) 
to  the  LOGGER. 

2.  L,  which  activates  the  LOGGER. 

3.  0,  which  deactivates  the  LOGGER. 

These  three  commands  allow  elementary  editing  capabilities  in  prepar¬ 
ation  of  pseudo-terminal  Inputs  for  batch  usage.  However,  an  on-site 
editor  will  be  required  to  produce  acceptable  batch  inputs  such  as 
Inserting  Calculate  connands,  changing  detection  criteria,  and  append¬ 
ing  several  LOGGER  productions  into  one  batch  Input  file. 

The  Interactive  commands  used  to  enter  and  exit  the  main  and  sub¬ 
menus  are  shown  In  Fig.  3. 
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VII.  ERRORS  AMD  DEFAULTS 

Default  values  are  always  displayed,  and  can  be  selected  by 
striking  the  CR  key. 

Four  types  of  errors  are  detectable  by  SNAPD. 

1.  Value  error,  e.g.,  a  negative  probability. 

2.  Syntax  error,  e.g.,  a  "7  ^  7"  numeric  string. 

3.  Input  error  indicated  by  a  nonzero  status  code. 

4.  I/O  format  errors  from  reading  incorrect  file  types. 

These  error  messages  always  refer  to  the  last  terminal  input. 

Some  errors  are  inherently  undetectable,  such  as  uninitialized 
wave  data.  The  user  must  take  great  care  to  avoid  them. 

Other  areas  where  the  user  must  exercise  caution  are  as  follows: 

1.  The  regional  grid  must  be  designed  to  break  at  ±180  deg  longi-* 
tude  and  be  indexed  starting  from  the  lower  left-hand  comer. 

2.  The  lowest  index  in  a  row  must  be  on  the  left-hand  border. 

3.  The  maximum  number  of  magnitude  corrections,  is  limited 

to  100. 

4.  For  the  master-event  and  S-P  calibration  sets,  and 
the  user  can  have  up  to  15  station  indices  per  epicenter. 

5.  Importance  sets  for  a  location  are  limited  to  a  maximtmi  of 
six. 
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VIII.  SAMPLE  RUN 


This  section  presents  an  abridged  printout  of  a  typical  Inter- 
active  session  for  a  SNAFO  run.  For  the  sake  of  brevity.  It  was 
necessary  to  omit  certain  portions  of  the  Input  and  output  printouts. 
Those  portions  of  the  printout  omitted  from  the  Interactive  Input 
session  are  typically  attenuation,  travel-time,  and  travel-time 
standard-deviation  table  entries  beyond  the  first  three  table  entries, 
and  station  Inputs  for  all  but  the  first  three  stations.  Detailed 
epicenter  output  Is  shown  for  only  the  first  three  epicenters. 
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Appendix 

DESCRIPTION  OF  THE  LOCATE  MODULE* 


GENERAL  DESCRIPTION 

The  module  LOCATE  predicts  the  error  a  seismic  network  would 
commit  In  locating  an  event  of  a  given  magnitude  and  true  location. 

The  expected  location  error  Is  expressed  In  the  form  of  confidence 
regions  for  the  epicenter  and  focal  depth  estimates  that  would  be 
derived  from  the  seismic  data  the  network  Is  predicted  to  observe. 

The  possible  data  consist  of  the  arrival  times  of  one  or  more  seismic 
phases  at  each  station  and  back  azimuth  measurements  obtained  from  Lg. 

The  location  confidence  regions  are  Inferred  from  the  following 
principal  Inputs: 

•  The  network  station  locations  with  respect  to  the  true  epi¬ 
center,  and  the  true  focal  depth  of  the  event. 

•  The  probability  that  each  seismic  phase  Is  observed  at  each 
station  (the 

e  The  standard  error  of  the  datum  obtained  from  each  detected 
phase  (the  or  • 

•  Travel  time  versus  distance  tables  contained  In  the  wave 
files  which  describe  the  theoretical  dependence  of  the  arrival¬ 
time  data  on  the  event  hypocenter. 

The  focal  depth  of  the  event  is  restricted  to  be  either  0,  13,  or  100  km, 
representing  an  explosion,  shallow  earthquake,  or  deep  earthquake,  re¬ 
spectively.  The  magnitude  of  the  event  Is  not  Input  explicitly,  but 
Is  reflected  In  the  observation  probabilities. 

3f 

Systems,  Science  and  Software  developed  the  LOCATE  module  and 
provided  this  general  description  of  Its  operation. 
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These  Inputs  determine  the  distribution  of  location  estimates 
that  would  be  derived  from  randomly  selected  sets  of  network  data; 
namely,  the  probability  density  function  f(£).  The  vector  £  denotes 
an  estimate  of  h,  the  four-element  vector  of  true  location  parameters: 

h  -  (8,  (>,  z,  t)^  , 

where  6  and  are  the  epicenter  coordinates,  z  Is  the  focal  depth, 
and  t  la  the  origin  time.  It  Is  assumed  that  h  is  the  "best  estl* 
mate  that  can  be  derived  from  any  sample  of  data;  e.g.,  the  minimum- 

A 

variance  unbiased  estimate.  The  distribution  of  h  can  be  shown  more 
completely  as  f(h|ti),  since  it  depends  on  the  true  location  h.  How- 

/K 

ever,  we  will  denote  It  simply  as  f(h)  for  the  sake  of  brevity. 

From  the  multivariate  distribution  f(h)  one  can  define  a  four- 

dlaenalonal  confidence  region  for  any  given  confidence  level  P^.  The 

location  estimate  derived  from  a  particular  network  data  set  would  have 

a  P  probability  of  falling  within  this  region, 
c 

For  nuclear  monitoring  applications  It  Is  useful  to  know  the 

/V  A 

separate  confidence  regions  of  the  epicenter  estimate  (9,  and  focal 
depth  estimate  z.  These  are  obtained  from  the  marginal  distributions 
f(9,  $)  and  f(z).  The  focal  depth  confidence  Interval,  z^  s:  z  ^  z^. 

Is  defined  by 


dz  f(z)  »  P^  , 


f(Zj^)  -  f(z2) 


The  two-dimensional  epicenter  confidence  region  E  and  Its  boundary 
3?  are  defined  by 
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// 


de  d<|)  fce. 


p 

c 


<^)  ~  constant  on  dE  . 


With  appropriate  asstanptlons  and  approximations,  f(h)  Is  a  normal 
distribution.  Implying  that  f(6,  (^)  and  f(z)  are  normal.  As  a  result, 
the  epicenter  confidence  region  Is  bounded  by  an  ellipse. 

The  principal  outputs  of  LOCATE  for  epicenter  j  are  A^ ,  the  area 
of  the  epicenter  confidence  ellipse,  and  d^ ,  the  length  of  the  focal 
depth  confidence  Interval.  In  addition,  LOCATE  returns  the  lengths 
and  orientations  of  the  principal  axes  of  the  epicenter  confidence 
ellipse. 

LOCATE  also  calculates  the  relative  Importances  of  the  Individual 
data  (l.e.,  station  1  observing  wave  k)  to  the  overall  network  loca¬ 
tion  capability.  The  Importance  of  a  datum  with  respect  to  a  location 
parameter  measures  how  much  the  confidence  region  for  that  parameter 
would  Increase  If  the  datum  were  removed  from  the  data  set  while  re¬ 
taining  all  the  other  data: 


Importance 


Confidence  region  without  datum 
Confidence  region  with  datum 


Two  Importances  are  defined  for  each  epicenter  J :  Che  importance 

of  station  1  and  wave  k  with  respect  to  the  epicenter  estimate,  and 

d  J 

'^Ik*  importance  with  respect  to  the  depth  estimate.  They  measure 

the  dependence  of  Aj  and  dj ,  respectively,  on  the  datum.  _ 

The  LOCATE  module  is  called  In  a  loop  In  which  the  true  event 
epicenter  Is  varied  from  call  to  call  while  the  true  focal  depth  and 
network  configuration  are  held  constant.  (The  detection  probabilities 
and  data  standard  errors  depend  on  the  event-station  distances,  so 
they  vary  with  the  epicenter.)  An  initialization  call  to  LOCATE  pre¬ 
cedes  the  loop,  and  a  finalization  call  follows  It. 
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Each  call  to  LOCATE  within  the  loop  returns  the  epicenter  and 
focal  depth  confidence  regions  for  the  event  having  the  current  epi¬ 
center.  The  data  importances,  however,  are  obtained  for  sets  of  epi¬ 
centers  r  ,  m  >  1,  M,  rather  than  for  each  individual  epicenter, 

n 

The  Importances  obtained  for  the  events  in  a  given  set  are  reduced 

to  six  values  for  each  statlon/wave  pair:  the  minimum,  maximum,  and 
A4  d-t 

average  values  of  and  over  the  set.  The  purpose  of  grouping 
the  importances  in  this  way  is  to  find  the  Importance  of  a  given 
station/wave  combination  over  user-specified  geographical  regions. 

DATA  MODELED 

LOCATE  models  data  from  five  seismic  phases:  P,  Pg,  Lg,  S,  and 
pP.  The  data  from  all  but  the  Lg  phase  are  assumed  to  be  arrival 
times;  Lg  is  assumed  to  provide  an  estimate  of  the  event-station  back 
azimuth.  The  P  wave  at  a  given  station  may  be  either  Pn,  mantle- 
refracted  P,  or  PKP,  depending  on  which  is  the  first  arrival  for  the 
event-station  distance.  Pg  is  considered  separately  since  in  practice 
it  can  be  observed  as  either  a  first  or  later  arrival.  The  S  wave  is 
also  a  composite  phase  (Sn,  S,  or  SRS) ,  but  it  can  be  used  only  from 
paths  for  which  the  P  to  S  velocity  ratio  (Vp/v.)  is  assumed  to  be 
calibrated  (i.e.,  the  S-P  technique). 

LOCATE 's  multiphase  capability  makes  it  possible  to  simulate  a 
variety  of  seismic  location  techniques  currently  in  use:  a  conven¬ 
tional  location  with  P  times,  the  Lg  back  azimuth  technique,  master- 
event  calibration  of  P  times,  origin-time  constraint  with  S-P  times 
from  Vp/Vg  calibrated  paths,  and  focal  depth  determination  with  pP-P 
times.  The  last  three  techniques  are  described  by  Evemden  [1969]. 
LOCATE  models  these  techniques  automatically  when  the  appropriate 
arrival  times  (P,  S,  or  pP)  belong  to  the  data  set  and  when  some  mod¬ 
ifications  are  made  to  the  arrival-time  standard  errors  as  discussed 
in  the  "ATSD"  subsection  of  Sec.  IV. 

LOCATE  calculates  the  confidence  regions  that  would  result  from 
a  simultaneous  estimation  of  all  four  location  parameters  (epicenter 
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coordinates,  focal  depth,  and  origin  time),  using  all  the  data  from 
the  various  phases  and  network  stations.  The  simultaneous  estimation 
approach  makes  most  efficient  use  of  the  data  and  produces  the  minimum- 
variance  estimate  of  each  location  parameter.  However,  this  approach 
is  not  always  taken  in  practice.  For  example,  the  S-P  and  pP-P  tech¬ 
niques  are  often  applied  separately  to  constrain  origin  time  and  focal 
depth,  and  then  the  Inferred  constraints  are  applied  in  an  epicenter 
determination  from  P  times.  This  and  other  variations  on  multiphase 
event  location  approximate  the  simultaneous  estimation  approach,  but 
they  are  not  explicitly  modeled  by  LOCATE. 

NETWORTH  [Wirth,  1977],  on  the  other  hand,  does  effectively  model 
the  step-wise,  as  opposed  to  simultaneous,  approach  to  epicenter  esti¬ 
mation.  NETWORTH  models  an  epicenter  determination  from  P  times,  given 
that  the  focal  depth  is  known  from  other  data  (e.g.,  constrained  by  S-P 
and  pP-P).  In  mathematical  terms,  NETWORTH  calculates  the  epicenter 
confidence  region  from  f(6,  *  z),  whereas  LOCATE  calculates  it  from 

f(6,  $)  which  is  marginal,  rather  than  conditional,  with  respect  to  z. 

STATISTICAL  FORMULATION  OF  EVENT  LOCATION 

To  properly  describe  the  location  capability  of  a  network,  the 
confidence  regions  determined  by  LOCATE  must  reflect  two  sources  of 
uncertainty  in  the  expected  network  data  set:  data  errors  and  un¬ 
detected  phases.  Finding  confidence  regions  that  account  for  the 
data  errors  is  a  straightforward  problem  in  linear  estimation.  Ac¬ 
counting  for  undetected  phases  is  a  more  difficult  problem. 

A  data  error  is  any  deviation  of  an  observed  datum  from  the  value 
that  would  be  predicted  theoretically  if  the  event  location  were 
known.  Thus,  the  observed  datum  from  the  V.th  seismic  phase  at  the 
i.th  station  is  related  to  the  true  event  location  h  by 

«ik  ‘  ^ik^^^  *  «ik  • 

where  g^^  is  the  datun,  (7^^^  is  the  "data  function"  expressing  the 
theoretical  dependence  of  the  datum  on  h,  and  e^j^  is  the  data  error. 
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Data  errors  can  result  from  observational  errors,  such  as  the  reading 
error  in  an  arrival-time  measurement,  or  modeling  errors,  such  as 
errors  in  the  travel-time  tables  used  for  predicting  theoretical 
arrival  times. 

With  some  simple  assumptions  about  e^  and  the  location 

confidence  regions  implied  by  a  given  set  of  detected  data  can  be 
easily  derived.  The  error  in  each  datun  is  assumed  to  be  a  normally 
distributed  zero-mean  random  variable  with  known  variance: 


-  0 


''"t'lk'  ■  °lk  • 


where  <7^^  is  the  standard  deviation  of  the  data  error  or,  more  simply, 
the  standard  error.  Then  if  each  is  approximated  by  a  linear 
function,  the  minimum-variance  unbiased  location  estimate  ti  is  a 
linear  function  of  the  observed  data.  Expressions  for  its  distribu- 

A 

tion  f(h)  and  its  confidence  taglons  are  rather  simple. 

Since  undetected  phases  contribute  no  data  to  the  network  data 
set,  the  subset  of  data  chat  would  be  available  for  deriving  a  loca¬ 
tion  estimate  is,  like  the  data  errors,  a  random  element  chat  must  be 
accounted  for  in  the  location  confidence  regions.  The  probability  of 
detection  p^  (suppressing  the  constant  index  j)  assigned  to  each  phase 
at  each  station  is  essentially  the  probability  that  the  (i,  V.)th 
datum  will  be  part  of  the  network  data  set.  Equivalently,  Che  (1,  k)th 
datum  has  a  finite  standard  error  with  probability  p^^  and  an  infinite 
standard  error  with  probability  (1  -  p^j^) . 

The  best  way  to  account  for  the  p^j^  is  not  obvious .  One  approach 
is  to  average  the  confidence  regions  computed  for  pseudo-randomly 
selected  subsets  of  data,  selected  according  to  Che  detection  prob¬ 
abilities.  Another  approach  is  to  compute  the  confidence  regions  for 
Che  most  probable  data  subset  containing  a  specified  number  of  data. 
These  approaches  are  implemented  In  METWORTH.  While  they  do  account 
for  the  effect  of  undetected  phases  on  location  performance,  they 
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define  the  expected  average  performance  of  a  network  in  a  rather  sub¬ 
jective  way.  Further,  the  first  approach  may  take  a  substantial 
amount  of  computation. 

The  most  meaningful  approach  to  the  problem  of  undetected  phases, 
and  the  approach  taken  by  LOCATE,  is  to  define  the  location  confidence 
regions  from  the  probability  distribution  of  h  that  treats  both  the 
data  subset  and  data  errors  as  random  variables.  That  is,  the  dis¬ 
tribution  of  location  estimates  that  would  result  from  random  sampling 
of  both  the  subset  of  detected  data  and  the  data  values  themselves. 

This  is  accomplished  by  treating  the  standard  errors  as  random  variables 
having  the  discrete  distribution 

Prob  -  p^^  , 

Prob  [o^  -  00]  .  1  -  p^^  , 

where  af,  is  the  standard  error  of  the  datum  obtained  from  a  detected 
ik 

phase.  The  distribution  described  earlier — namely,  that  of  the 
minimum-variance  estimate  obtained  from  a  given  data  subset — is  ac¬ 
tually  the  oonditionaZ  distribution  f(£|<7^j^).  The  distribution  that 
accounts  for  undetected  phases  is  rrarginoL  with  respect  to  the 
We  now  show  the  latter  as  f (h) . 

The  exact  marginal  distribution  of  ti  is  quite  difficult  and 
costly  to  compute.  Further,  it  is  not  a  normal  distribution,  even 
when  the  errors  are  normal,  Aince  the  standard  errors  are  not  normally 
distributed.  Therefore,  obtaining  the  none Hips oldal  confidence  re- 
glons  from  f(h)  would  be  quite  involved  even  if  the  distribution  were 
known. 

A 

Fortunately,  an  approximation  to  the  marginal  distribution  of  h 
can  be  evaluated  as  easily  as  its  conditional  distribution.  The 
approximation  in  fact  is  to  evaluate  the  conditional  distribution  with 
weighted  standard  errors: 
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f(h)  -  f(h|o^^  -  , 


where 


Thus,  the  effect  of  undetected  phases  can  be  approximated  by  assuming 
every  phase  Is  detected  and  produces  a  datum  whose  error  Is  Increased 
by  a  factor  depending  on  the  probability  of  detection. 

The  approximate  f(h)  Is  a  normal  distribution,  so  the  approxi¬ 
mate  confidence  region  of  ti  Is  ellipsoidal.  The  accuracy  of  the 
approximation.  It  turns  out.  Is  Indicated  by  the  size  of  the  data 

importances.  The  approximation  is  accurate  when  the  Individual  data 

Aj  d4 

are  not  very  important  and  >  1)  • 
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